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ABSTRACT  OF  THE  DISSERTATION 


Velocity  and  Attenuation  Profiles  in  the  Monterey  Deep-Sea  Fan 


by 


Richard  K.  Brienzo 


Doctor  of  Philosophy  in  Electrical  Engineering  (Applied  Ocean  Science) 
University  of  California,  San  Diego,  1987 
Professor  William  S.  Hodgkiss,  Chair 


Data  obtained  during  a  refraction  experiment  is  used  to  estimate  velo¬ 
city  and  attenuation  profiles  in  an  area  of  thick  sediments  (2.5  -  3  km)  on  the 
Monterey  Deep-Sea  F an.  A  20  element  vertical  hydrophone  array  was  deployed 
at  mid-depth  in  2800  meters  of  water.  Explosive  sources,  set  for  a  depth  of  1820 
meters,  were  detonated  at  ranges  between  3.5  km  to  37  kilometers  from  the 
array.  Estimates  of  velocity  as  a  function  of  depth,  and  attenuation  as  a  func¬ 
tion  of  frequency  and  depth  are  obtained  from  an  analysis  of  the  pressure  time 
series  generated  by  the  explosive  charges  and  received  at  the  array. 

To  find  the  velocity  profile,  the  sediment  is  modeled  as  a  horizontally 
layered,  laterally  homogeneous  medium.  A  least  squares  solution  is  found  for 
the  velocity  gradients  in  each  layer  of  the  model.  Velocity  as  a  function  of 
depth  is  obtained  by  integrating  these  gradients.  A  second  approach  to  inferring 
velocity  structure  utilizes  linear  programming.  Errors  in  the  data  are  easily 
incorporated  into  the  problem  using  this  formulation.  The  method  takes  upper 
and  lower  bounds  on  the  input  data  and  gives  as  a  solution  upper  and  lower 
bounds  on  the  velocity  profile.  All  velocity  profiles  that  are  consistent  with  the 
data  lie  within  these  bounds. 

A  method  of  spectral  ratios  is  used  to  estimate  attenuation  in  the  sedi¬ 
ment  as  a  function  of  frequency  and  depth.  The  sediment  is  again  modeled  as  a 
layered  medium,  with  each  layer  having  a  separate  attenuation  coefficient. 
Solving  for  each  of  the  coefficients  gives  attenuation  as  a  function  of  depth. 
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1.  Introduction 


Both  the  water  and  bottom  sediment  act  as  transmission  media  for 
acoustic  signals.  Since  there  is  a  general  increase  in  sound  velocity  with  depth 
in  the  sediment,  much  of  the  acoustic  energy  that  enters  the  sediment  is 
refracted  back  into  the  water  column.  Knowing  the  velocity  and  attenuation  of 
sound  in  the  sediment  is  therefore  important  in  modeling  sound  propagation  in 
the  ocean  as  well  as  providing  information  about  sediment  type  and  structure. 
In  this  report,  propagation  of  low  frequency  sound  in  the  sediments  will  be  stu¬ 
died  through  the  estimation  of  velocity  and  attenuation  profiles  in  the  sediment. 

In  April  1981,  the  Marine  Physical  Laboratory  conducted  a  refraction 
experiment  in  an  area  of  thick  sediments  (2.5-3  km)  on  the  Monterey  fan.  The 
research  vessel  FLIP  (FLoating  Instrument  Platform)  was  moored  in  2800  meters 
of  water  and  utilized  as  a  receiving  station.  A  20  element  vertical  hydrophone 
array  with  equal  sensor  spacing  and  a  total  aperture  of  475  meters  (25  meter 
spacing  between  elements)  was  deployed  at  mid-depth  in  the  water  column. 
Explosive  sources,  set  for  a  depth  of  1820  meters,  were  detonated  at  ranges 
between  3.5  km  to  37  km  from  FLIP. 

Acoustic  energy  travels  along  numerous  paths  in  going  from  the  shot  to 
the  hydrophone  array.  Of  interest  in  this  study  is  acoustic  energy  which  enters, 
propagates  through  the  sediment,  and  is  subsequently  refracted  back  into  the 
water  column  before  being  received  at  the  array.  Both  velocity  as  a  function  of 
depth,  and  attenuation  as  a  function  of  frequency  and  depth  are  obtained  from 
an  analysis  of  the  pressure  time  series  generated  by  the  explosive  charges  and 
received  at  the  array. 

The  report  is  organized  as  follows.  Chapter  2  contains  a  description  of 
the  refraction  experiment  on  the  Monterey  Fan.  Examples  of  the  array  time 
series  are  given  and  individual  arrivals  in  the  time  series  are  related  to  specific 
raypaths.  The  chapter  concludes  with  a  discussion  of  some  of  the  characteristics 
of  shot  spectra. 

In  Chapter  3,  the  array  time  series  are  used  to  estimate  velocity  in  the 
sediment  as  a  function  of  depth.  Two  methods  are  used  to  obtain  information 
about  the  velocity  structure.  The  sediment  is  modeled  as  a  horizontally  layered, 
laterally  homogeneous  medium.  A  least  squares  solution  is  found  for  the  velocity 
gradients  in  each  layer  of  the  model.  Velocity  as  a  function  of  depth  is  obtained 
by  integrating  these  gradients.  Inversion  for  the  velocity  profile  by  this  method 
is  discussed  in  the  context  of  a  general  linear  inverse  problem  in  Appendix  Al.  A 
second  approach  to  inferring  velocity  structure  utilizes  linear  programming. 
Errors  in  the  data  are  easily  incorporated  into  the  problem  using  this  formula¬ 
tion.  The  method  takes  upper  and  lower  bounds  on  the  input  data  and  gives  as 
a  solution  upper  and  lower  bounds  on  the  velocity  profile.  All  velocity  profiles 
that  are  consistent  with  the  data  lie  within  these  bounds. 

A  method  of  spectral  ratios  is  used  in  Chapter  4  to  estimate  attenua¬ 
tion  in  the  sediment  as  a  function  of  frequency  and  depth.  Spectral  ratios  are 
formed  by  dividing  the  spectrum  of  an  arrival  which  travels  through  the 


1 


R.  K.  Brienzo 


sediment  by  that  of  an  arrival  which  travels  only  through  the  water  column.  A 
derivation  in  Appendix  A2  shows  that  attenuation  in  the  sediment  is  propor¬ 
tional  to  the  slope  of  the  spectral  ratio.  The  sediment  is  again  modeled  as  a 
layered  medium,  with  each  layer  having  a  separate  attenuation  coefficient. 
Solving  for  each  of  the  attenuation  coefficients  gives  attenuation  as  a  function  of 
depth. 
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2.  Data 


2.1  Experiment 

In  April  1981,  the  Alnrine  Physical  Lnborntory  conducted  a,  refrsiction 
experiment  in  an  area  of  thick  sediments  (2.5-3  km)  on  the  Monterey  fan 
{ZG°  37.31  '.V,  122"  41.76' VT)  .  The  research  vessel  FLIP  (FLoating  Instrument 
Platform)  was  moored  in  2800  meters  of  water  and  utilized  as  a  receiving  sta¬ 
tion  (Figure  2  1).  Figure  2.2  shows  the  general  area  of  the  experiment. 

A  20  element  vertical  hydrophone  array  with  equal  sensor  spacing  and 
a  total  aperture  of  475  meters  (25  meter  spacing  between  elements)  was 
deployed  at  mid-depth  in  the  water  column.  The  hydrophones  were  numbered 
from  1  to  20,  with  number  1  being  the  bottom  element.  Hydrophone  number  1 
was  located  at  1912  meters,  which  placed  the  top  element  at  1437  meters,  and 
the  midpoint  of  the  array  at  1675  meters.  Hydrophones  11-20  were  set’ at  a 
gain  of  approximately  40  dB  higher  than  1  -  10  to  accommodate  a  wide 
dynamic  range  of  signals.  Hydrophones  1  and  4  were  later  found  to  have  mal¬ 
functioned  and  are  not  used  in  the  data  analysis.  The  analog  hydrophone  out¬ 
puts  were  lowpass  filtered  at  400  Hz,  digitized  at  1  kHz  per  channel  and 
recorded  on  magnetic  tape. 

Explosive  sources  were  launched  over  two  separate  seven  hour  periods 
at  ranges  from  3.5  km  to  37  km  from  FLIP.  Figure  2.3  indicates  the  tracks  of 
the  two  shot  runs.  Both  shallow  (300'  Mk  82  SUS)  and  deep  (6000'  Mk  94  SUS  + 
256#TNT)  charges  were  utilized.  Only  the  deep  charges  were  used  for  the  work 
presented  in  this  report.  In  addition  to  generating  more  energy,  they  are  far 
enough  from  bottom  and  surface  reflecting  layers  so  that  reflected  energy  does 
not  interfere  with  energy  arriving  along  other  paths  of  interest. 

The  position  of  the  shooting  ship  was  monitored  from  FLIP  by  radar. 
Water  temperature  data  were  collected  using  1800  meter  XBT’s  from  the  shoot¬ 
ing  ship  during  both  shot  runs.  A  total  of  19  temperature  profiles  were  collected 
over  a  period  of  two  days. 

The  Monterey  fan  (Figure  2.2)  has  been  the  subject  of  many  studies  [27, 
47,  48,  49,  67].  It  is  primarily  fed  by  the  Monterey  and  Ascension  canyon  sys¬ 
tems.  The  sediments  are  mainly  terrigenous,  consisting  of  very  fine  sands  in  the 
valleys,  and  mostly  silts  and  clays  away  from  the  valleys  and  channels.  The  sea 
floor  in  the  area  of  this  experiment  is  relatively  flat  (Figure  2.3).  The  depth  was 
2800  meters  at  FLIP,  increased  to  a  depth  of  2930  meters  at  a  range  of  20.2 
kilometers,  and  to  3230  meters  depth  at  a  range  of  32.4  kilometers. 

Figure  2.4  illustrates  a  few  of  the  raypaths  that  sound  may  take  in 
traveling  from  the  shot  to  the  array.  Paths  that  interact  with  the  bottom  more 
than  once  will  not  be  considered  in  this  analysis.  This  figure  also  shows  the  cov- 
erage  (sampling)  of  the  sediment  by  each  type  of  ray.  Bottom  interacting  rays 
with  a  small  angle  of  incidence  (BR)  graze  the  sediment  while  the  steeper  ray- 
paths  (S-BR-S)  sample  deeper  portions  of  the  sediment.  The  maximum  depth  in 
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Figure  2.1  Refraction  experiment  on  the  Monterey  Fan. 

A  20  element  hydrophone  array  was  deployed  at  mid-depth  in  2800  meters  of  water.  . 
Explosive  charges  were  detonated  at  ranges  between  3.5  to  37  km  from  FLIP. 
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Figure  2.2  The  Monterey  Deep  Sea  Fan. 

Depths  are  reported  in  meters  x  10 
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Figure  2.3  Location  of  the  refraction  experiment. 

Depths  are  reported  in  meters  x  10'^. 
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Figure  2.4  Raypaths  from  shot  to  array. 
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the  sediment  for  which  arrivals  in  this  experiment  are  still  discernible  is  about 
625  meters. 


Ravpaths  (see  Fteure  2.4) 

label 

path 

D 

direct  water  path 

S 

surface  reflection 

B 

bottom  reflection 

BR 

bottom  refraction 

BR-S 

bottom  refraction  -  surface  reflection 

S-BR 

surface  reflection  -  bottom  refraction 

S-BR-S 

...  surface  reflection  -  bottom  refraction  -  surface  reflection 

2.2  Time  Series 

Figures  2.5a,  2.5b,  and  2.5c  are  examples  of  the  type  of  data  that  were 
received  at  the  hydrophone  array.  Each  figure  shows  the  output  of  hydrophones 
f  1  ■  20  for  a  single  shot.  The  2.048  seconds  of  data  plotted  here  contains  several 
distinct  arrivals.  Each  of  the  arrivals  represent  acoustic  energy  that  has  trav¬ 
eled  along  a  different  path.  The  arrivals  in  these  figures  have  been  labeled  to 
correspond  to  the  raypaths  in  Figure  2.4.  The  method  used  to  associate  an 
arrival  at  the  array  with  a  raypath  is  described  in  Section  3.3.1. 

The  data  presented  in  these  three  figures  are  from  the  high  gain 
(upper)  half  of  the  array.  The  larger  amplitude  arrivals  have  saturated,  and  no 
inference  should  be  made  about  the  relative  amplitudes  of  arrivals.  These 
figures  are  meant  to  show  angles  at  which  the  arrivals  come  into  the  array,  and 
the  time  relationships  between  them.  Acoustic  energy  coming  from  above  the 
array  arrives  at  the  top  hydrophone  (20)  first,  and  the  line  of  arrivals  appears  to 
be  slanted  to  the  left  (eg.  shot  43,  BR-S).  Energy  coming  from  the  bottom 
arrives  at  the  bottom  element  first,  and  the  arrivals  appear  to  be  slanted  to  the 
right  (eg.  shot  43,  S-BR).  At  short  ranges,  the  arrivals  are  all  separated  in  time. 
As  the  shot  range  increases,  differences  in  the  path  lengths  (hence  differences  in 
time  between  the  arrivals)  decrease,  and  the  arrivals  begin  to  run  together. 

Marine  refraction  experiments  commonly  have  the  source  and/or 
receiver  near  the  surface  or  bottom.  This  causes  sediment  refracted  arrivals  to 
be  contaminated  by  reflections  from  these  interfaces.  Refraction  experiments 
where  the  source  is  near  the  surface  cannot  be  carried  out  in  deep  water  since 
the  simple  raypaths  refract  completely  in  the  water  column  and  never  enter  the 
sediment.  In  the  Monterey  fan  experiment  both  source  and  receiver  were  placed 
at  mid-depth  in  relatively  deep  water.  At  short  ranges  (Fig.  2.5a),  the  arrivals 
are  well  separated  in  time.  As  the  shot  range  increases,  arrivals  begin  to  run 
together  (Fig.  2.5c).  By  18  km,  the  bottom  refraction  (BR)  disappears  and  the 
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Figure  2.5a  Array  Time  Series  -  Shot  15 

Arrivals  marked  D,  B  and  S  correspond  to  the  raypaths  of  Figure  2.4. 
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Figure  2.5b  Array  Time  Series  -  Shot  43 

Arrivals  marked  D,  BR,  S,  BR-S  and  S-BR  correspond  to  the  raypaths  of  Figure  2.4. 
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Figure  2.5c  Array  Time  Series  -  Shot  25 

Arrivals  marked  D,  BR,  S,  BR-S,  S-BR  and  S-BR-S  correspond  to  the 
raypaths  of  Figure  2.4. 
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ray  refracts  in  the  water  column.  The  other  sediment  refracted  arrivals  (BR-S, 
S-BR,  S-BR-S)  exist  and  are  separable  out  to  about  25  km. 

The  time  series  in  Figure  2.5  also  indicate  that  bottom  refracted  energy 
does  not  appear  at  short  ranges.  When  a  wavefront  encounters  a  change  in 
acoustic  impedance,  part  of  the  energy  will  be  reflected  and  the  remaining  part 
transmitted  At  short  ranges,  the  angle  of  incidence  at  the  array  is  large  and 
the  refracted  energy  penetrates  deeply  in  the  sediment.  Rays  with  very  deep 
turning  points  have  long  path  lengths  and  are  so  highly  attenuated  that  they  do 
not  appear  in  the  time  series.  As  the  range  increases,  the  angle  of  incidence  at 
the  array  decreases,  resulting  in  shallower  turning  points.  The  path  length 
becomes  shorter,  decreasing  the  amount  of  attenuation,  and  refracted  energy 
begins  to  appear.  Figure  2.6  (reprinted  from  reference  [7])  shows  the  rapid 
decrease  in  turning  point  depth  with  increasing  range.  Reflections  off  the  sea 
floor  lose  much  of  their  energy  and  appear  as  low  amplitude,  high  frequency 
arrivals. 

Figure  2.7  is  a  plot  of  the  "stack",  and  illustrates  some  of  the  above 
features.  Each  time  series  in  this  figure  is  the  output  of  a  single  channel 
(number  15)  from  a  different  shot.  Again,  note  that  the  larger  amplitude 
arrivals  have  saturated,  so  the  time  series  presented  here  contain  no  information 
about  relative  amplitudes.  The  bottom  time  series  is  from  a  shot  with  a  range  of 
5741  meters.  There  are  only  three  distinct  arrivals:  a  direct  water  path  (D),  a 
bottom  reflection  (B),  and  the  surface  reflection  (S).  By  the  time  the  range  has 
increased  to  7871  meters,  a  bottom  refraction  (BR)  has  emerged.  At  a  range  of 
9075  meters,  the  first  BR-S  arrival  is  seen.  Similarly,  the  S-BR-S  arrival  does 
not  become  apparent  until  a  range  of  16,112  meters.  The  deepest  ray  seen  in 
this  data  has  a  turning  depth  of  about  625  meters. 


2.3  Spectra 

Explosive  sources  generate  energy  over  a  wide  range  of  frequencies. 
Upon  detonation  the  explosive  material  is  rapidly  converted  to  a  gas,  and  a 
shock  wave  with  a  very  steep  front  is  radiated  into  the  water  If  a  pure  impulse 
was  generated,  an  explosive  source  would  have  a  flat  spectrum  This  is  not  quite 
the  case,  however.  As  the  gas  bubble  expands,  internal  pressure  decreases,  but 
because  of  inertia,  the  gas  bubble  overshoots  its  equilibrium  radius.  The  pres¬ 
sure  in  the  gas  bubble  goes  slightly  negative,  and  the  bubble  begins  to  contract. 
The  gas  bubble  again  overshoots  its  equilibrium  radius,  the  pressure  becomes 
positive  and  the  bubble  begins  to  expand  again.  The  amplitude  of  these  bubble 
pulses  decays  fairly  quickly,  and  no  more  than  four  or  five  pulses  are  usually 
seen.  The  decaying,  oscillatory  nature  of  the  resulting  waveform  produces  a 
spectrum  that  is  broadband,  but  not  flat.  Characteristics  of  the  shot  spectrum 
will  be  discussed  below.  Various  relationships  between  the  bubble  pulse  period, 
charge  weight,  shot  depth,  and  other  parameters  may  be  found  in  [5.  381. 
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Figure  2.6  Refracted  paths  as  a  function  of  range 

As  range  increases,  the  depth  of  the  turning  point  decreases.  The  arrow  indicates  the  ray  where 
refracted  energy  is  first  seen  at  the  array.  (Reprinted  from  Christensen,  R.E.,  Frank,  J.A.,  and 
Geddes,  W.H.,  "Low  frequency  propagation  via  shallow  refracted  paths  through  deep  ocean 
unconsolidated  sediments",  JASA,  vol.57,  No.6,  1421-1426,  June  1975] 
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Figure  2.7  Received  time  series  as  a  function  of  range 

Each  time  series  is  a  single  channel  from  a  different  shot.  Shot  range 
increases  from  bottom  to  top.  Note  that  the  ranges  are  not  equally  spaced. 
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Figure  2.8  shows  the  output  of  channel  3  from  shot  41.  Two  refracted 
arrivals  are  present.  There  is  a  relatively  large  BR  arrival,  and  a  very  deep  and 
thus  highly  attenuated  BR-S  arrival.  The  individual  arrivals  plotted  below  are 
scaled  to  have  the  same  maximum  amplitudes.  Spectra  for  these  four  arrivals 
are  shown  in  the  top  of  Figure  2.9.  These  were  found  by  separating  the  indivi¬ 
dual  arrivals  of  each  type  (D,  BR,  S,'BR-S)  from  the  eight  good  channels  of  the 
low  gam  half  of  the  array.  Discrete  Fourier  transforms  of  each  type  were  calcu¬ 
lated  and  incoherently  averaged.  The  y-axis  in  these  plots  is  20  log  (magni¬ 
tude).  The  time  series  was  not  calibrated,  so  in  this  case  amplitudes  represent 
sampled  A-D  levels.  Calibration  is  not  necessary  since  it  is  ratios  of  the  various 
spectra  that  are  of  interest. 

Spectra  for  the  D,  BR  and  BR-S  paths  are  compared  in  the  bottom  plot 
of  Figure  2.9.  The  direct  path  arrival  (D)  has  a  relatively  broadband  spectrum 
with  somewhat  regularly  spaced  peaks.  When  compared  to  the  direct  path  spec¬ 
trum,  the  bottom  refracted  (BR)  spectrum  appears  to  fall  off  more  rapidly  with 
increasing  frequency,  being  nearly  equal  to  the  direct  path  spectrum  at  the  bub¬ 
ble  pulse  frequency  (50  Hz),  and  decreasing  by  about  8  dB  at  350  Hz.  The  BR-S 
spectrum  also  shows  an  attenuation  of  the  higher  frequencies.  This  decrease 
with  respect  to  the  direct  path  spectrum  will  be  exploited  to  find  attenuation  as 
a  function  of  frequency  and  depth  in  Chapter  Four 

A  simple  way  to  view  the  gross  features  of  a  shot  spectrum  is  to  con¬ 
sider  the  arrival  as  being  a  series  of  impulses.  If  it  were  an  impulse  train  with  a 
time  spacing  of  T  seconds,  then  its  Fourier  transform  would  be  an  impulse  train 
with  a  frequency  spacing  of  l/T.  Although  a  shot  arrival  does  not  consist  of  a 
periodic  impulse  train  -  it  is  of  finite  length,  has  decreasing  amplitude,  and  the 
bubble  pulse  period  is  not  constant  -  there  is  a  rough  correspondence  between 
the  direct  arrival’s  primary  bubble  pulse  period  of  20  msec,  and  a  spacing  of 
about  50  Hz  in  the  peaks  of  its  spectrum 

This  is  illustrated  in  Figure  2.10,  where  the  direct  path  arrival  is 
represented  by  four  impulses  with  spacing  and  amplitude  matching  the  peaks  in 
the  time  series.  The  spectrum  of  this  impulsive  model  along  with  the  spectrum 
of  the  direct  path  arrival  is  plotted  in  Figure  2.10  (bottom).  Although  it  is  a 
simple  model,  its  spectrum  looks  very  similar  to  that  of  the  actual  arrival. 

More  complicated  ways  to  model  the  shot  spectrum  have  been  devised. 
These  models  generally  make  the  impulses  look  more  like  the  real  time  series  - 
for  instance,  putting  increasing  and  decreasing  exponentials  on  either  side  of  the 
impulses.  Modeling  the  shot  spectrum  is  not  a  concern  in  this  work.  There  is 
little  attenuation  of  low  frequency  sound  traveling  short  distances  through  the 
water  column  [36,  65  (pg.  106)),  and  the  direct  path  arrival  (D)  is  a  good  approx¬ 
imation  of  the  source  spectrum. 

The  array  time  series  provides  the  basic  data  set.  In  Chapter  Three, 
sediment  refracted  arrivals  from  this  time  series  will  be  used  to  find  velocity  as 
a  function  of  depth  in  the  upper  sediments.  In  Chapter  Four,  differences 
between  the  direct  path  (source)  and  sediment  refracted  spectra  will  be  used  to 
infer  attenuation  as  a  function  of  frequency  and  depth  in  the  sediment. 
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Figure  2.8  Single  channel  time  series. 

Individual  arrivals  plotted  below  are  scaled  to  have  the  same  maximum  amplitude. 
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nuiica  refUcsd  rprcouiD 


bottom  refracted  -  surface  reflected  T-rciJj-j=i 


Figure  2.9  Spectra  of  the  arrivals  from  Figure  2.8. 
y-axis  is  20*log(magnitude) 
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A  T1  =  20  msec 
A  T2  =  14  msec 
A  T3  =  12  msec 


A1  =  1. 
A2  =  .82 
A3  =  .29 
A4  =  .1 


spectrum  of  the  model 
spectrum  of  the  shot  arrival 


Figure  2.10  Impulsive  model  of  the  shot  spectrum. 

Top  figure  is  a  direct  path  arrival  (D).  Middle  figure  is  an  impulsive  model. 
Bottom  figure  shows  the  spectra  of  both.  The  spectrum  of  the  model  has  been 
slightly  offset  to  separate  the  two  curves. 
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3.  Velocity  Profile  Estimation 


3.1  Introduction 

In  this  chapter  an  estimate  of  the  velocity  profile  in  the  upper  sediment 
(to  about  625  meters  depth)  is  found.  Since  much  of  the  acoustic  energy  enter¬ 
ing  the  sea  floor  is  refracted  back  into  the  water  column,  the  velocity  profile  is 
important  for  modeling  sound  propagation  in  the  ocean.  An  accurate  velocity 
profile  will  also  be  required  in  Chapter  Four  for  use  in  estimating  attenuation  as 
a  function  of  depth  in  the  sediment. 

An  overview  of  what  is  known  about  the  velocity  structure  in  the  upper 
sediments  is  first  presented.  Next,  steps  taken  to  convert  the  received  time  series 
at  the  hydrophone  array  into  sediment  traveltimes  (horizontal  range  in  the  sedi¬ 
ment  as  a  function  of  time)  are  described.  This  procedure  involves  several  steps. 
A  raypath  must  first  be  associated  with  each  arrival  at  the  array.  Shot  range 
and  depth,  and  depth  of  the  sea  floor  must  be  estimated.  Finally,  the  range  and 
time  spent  by  the  ray  in  the  water  column  must  be  subtracted  from  the  total 
range  and  time,  yielding  range  and  time  of  the  ray  in  the  sediment.  This  set  of 
ranges  and  times  form  the  traveltime  data,  giving  time  as  a  function  of  range 
T(X).  After  re-parameterizing  the  data  into  functions  of  ray  parameter  p,  they 
are  then  inverted  to  yield  velocity  as  a  function  of  depth. 

Two  methods  are  used  to  obtain  information  about  velocity  as  a  func¬ 
tion  of  depth.  A  least  squares  solution  will  be  used  to  find  velocity  gradients  in 
a  layered  model.  A  velocity  profile  is  obtained  from  these  gradients.  A  second 
approach  will  pose  the  problem  as  a  linear  programming  problem.  Errors  in  the 
data  are  easily  incorporated  into  the  problem  using  this  formulation.  This 
method  takes  upper  and  lower  bounds  on  the  input  data  and  gives  as  a  solution, 
upper  and  lower  bounds  on  the  velocity  profile.  All  velocity  profiles  that  are 
consistent  with  the  data  lie  within  these  bounds. 

Placing  meaningful  bounds  on  the  solution  is  important  in  interpreting 
the  result.  The  least  squares  solution  provides  an  answer;  it  gives  a  velocity 
profile  which  can  then  be  used  for  modeling  sound  propagation.  Because  of 
finiteness  of  the  data  set,  and  errors,  there  are  an  infinite  number  of  models 
which  fit  the  data.  Unless  a  statement  regarding  the  quality  of  the  result  can  be 
made,  all  that  can  be  said  is  that  the  solution  fits  the  data  while  minimizing  a 
particular  error  criterion  (the  /2  norm).  Although  rather  pessimistic  bounds  are 
produced  by  the  linear  programming  method,  all  models  that  are  consistent  with 
the  data  lie  within  these  bounds.  These  bounds  may  be  used  to  indicate  the  reli¬ 
ability  of  the  least  squares  solution  and  show  what  can  be  inferred  about  the 
velocity  structure  from  a  finite  data  set  containing  errors.  Using  these  two 
different  approaches  will  provide  a  clearer  interpretation  of  the  results. 
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3.2  Velocity  structure  in  the  upper  sediments 

Hamilton,  in  a  series  of  papers  [16,  21,  22,  24,  25,  26],  has  summarized 
much  of  what  is  known  about  the  velocity  structure  in  sediments.  In  general, 
the  velocity  of  sound  in  sediments  increases  with  increasing  depth.  Velocity  gra¬ 
dients  are  usually  highest  near  the  surface,  decrease  with  depth,  and  approach  a 
constant  by  about  700  meters.  Mean  values  for  velocity  gradients  for  all  sedi¬ 
ments  vary  from  1.28  sec at  the  surface  to  .6  sec “*  below  700  meters. 


8.2.1  Causes  of  velocity  gradients 


Reduction  of  porosity 

The  most  important  factor  is  a  reduction  of  porosity  with  increasing 
depth  in  the  sediment.  Porosity  is  the  volume  of  water-filled  pore  space  per 
volume  of  sediment.  The  amount  of  pore  space  in  sediment  is  related  to  the 
size,  shape,  distribution,  mineralogy  and  packing  of  the  solid  grains.  In  general, 
porosity  increases  as  grain  size  decreases.  As  small  particles,  such  as  silts  and 
clays  fall  and  touch  each  other,  they  are  loosely  held  together  by  interparticle 
forces  and  form  a  three-dimensional  "cardhouse"  structure.  These  sediments 
(for  example,  turbidites,  clays  and  silts)  near  the  surface  are  highly  porous,  typi¬ 
cally  85%.  Sands,  which  are  larger  and  heavier,  have  porosities  between  35  to 
50%.  As  overburden  pressure  increases  (with  increasing  depth),  sediment  parti¬ 
cles  are  forced  closer  together,  decreasing  the  porosity.  As  the  number  of  inter- 
particle  contacts  increase,  sliding  friction  between  particles  increases.  This 
causes  an  increase  in  rigidity,  which  leads  to  an  increase  in  velocity.  Porosity 
reduction  accounts  for  about  66%  of  the  total  gradient. 


Increase  in  temperature 

An  increase  in  temperature  with  depth  is  caused  by  heat  flow  through 
the  sediment  from  the  earth  s  crust  and  mantle.  Hamilton  i22  cites  a  change  of 
about  25  *  G  as  depth  in  the  sediment  increases  from  0  to  500  meters.  Sound 
velocities  vary  with  temperature  in  sediments  about  the  same  way  as  they  do  in 
sea  water;  the  velocity  of  sound  increases  with  increasing  temperature.  This 
accounts  for  about  17%  of  the  total  gradient. 


Increases  in  mineral  frame  rigidity  due  to  lithification 

As  sediment  particles  are  forced  together  under  increasing  pressure, 
they  come  into  closer  contact.  This  aids  in  cementation  of  the  particles 
(lithification)  and  increases  the  rigidity  of  the  mineral  frame,  which  increases  the 
velocity.  At  depths  below  about  300—600  meters  the  sediment  is  semilithified  or 
lithified  and  is  called  a  mudstone  or  claystone  After  lithification  there  is  little 
decrease  in  porosity  with  increasing  depth.  Lithification  is  responsible  for  about 
15%  of  the  average  gradient. 
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It  should  be  noted  that  there  is  not  a  smooth  transition  between 
unlithified  and  lithified  sediments.  For  example,  calcareous  sediments  easily 
lithify  into  relatively  high  velocity  chalk  and  limestone.  These  may  appear  at 
shallow  depths  (150—300  meters),  but  they  are  often  interbedded  with  soft  lower 
velocity  unlithified  sediments.  This  can  cause  scatter  around  the  velocity-depth 
function. 

For  terrigenous  sediments,  the  average  linear  velocity  gradient  between 
the  surface  and  500  meters  depth  is  about  1  sec-1.  Gradients  may  be 
significantly  lower  in  areas  where  sediments  are  rapidly  deposited  than  in 
regions  farther  from  the  source.  For  example,  Hamilton  [22,  pg.  920]  reports  a 
gradient  of  .86  near  the  source  of  the  Bengal  Fan,  and  a  value  of  1.87  farther  to 
the  South.  When  sediments  accumulate  in  thick  sections  the  porosity  does  not 
have  time  to  be  fully  reduced  under  overburden  pressure.  Since  2/3  of  the  aver¬ 
age  velocity  gradient  in  silt-clays  and  turbidites  is  caused  by  a  decrease  of 
porosity  with  increasing  depth,  areas  where  the  sedimentation  rate  is  high  have 
the  smallest  velocity  gradients.  Hamilton  [24,  pg.  1354]  reports  upper  sediments 
(0-100  meters)  for  5  shelf  basins  to  have  gradients  between  .3  to  .9  sec-1,  with  an 
average  gradient  of  .7  sec-1. 


3.2.2  Velocity  anisotropy 

Media  having  velocities  parallel  to  layers  that  are  different  from  veloci¬ 
ties  perpendicular  to  layers  are  called  anisotropic.  If  the  cardhouse  structure 
were  to  collapse  (under  sufficient  pressure)  to  form  a  parallel-oriented  section, 
velocity  anisotropy  may  be  induced.  Hamilton  claims  there  should  be  little 
anisotropy  from  this  effect  to  about  200  meters,  and  a  5  to  10%  increase  in  velo¬ 
city  parallel  to  the  sea  floor  at  400  to  600  meters  depth  [16,  pg.  4442-44431. 
Another  cause  of  anisotropy  that  may  be  more  significant  in  regions  of  high 
deposition  (such  as  deep-sea  fans),  is  alternating  parallel  layers  of  different 
materials. 


3.2.3  Velocity  at  the  water-sediment  interface 

Several  authors  have  indicated  that  the  velocity  in  the  upper  layer  of 
some  sediments  is  slightly  lower  than  that  in  the  water  column.  Hamilton 
reports  a  ratio  of  sediment  to  bottom  water  velocity  of  .984  for  pelagic  clay  and 
a  value  of  .993  for  turbidites  [22,  pg.  920],  The  turbidite  samples  included  silts 
and  clays. 

The  velocity  ratio  at  36"  31'  N,  123"  17'  W  [12]  was  reported  to  be  .988. 
This  is  very  close  to  the  site  of  the  Monterey  fan  refraction  experiment  (which 
was  at  36"  37.31'  N,  122"  41.76'  W)  and  the  value  of  .988  will  be  used  in  the  sub¬ 
sequent  analysis.  The  mean  depth  of  the  sea  floor  will  be  taken  to  be  2853 
meters.  The  velocity  of  sound  in  the  water  at  2853  meters  is  1502.8  meters/sec, 
so  the  velocity  of  sound  in  the  upper  layer  of  sediment  (interface)  is 

(1502. 8meters /second)  X  .988  =  1484.8 meters / second. 

Velocities  in  the  water  column  are  based  on  results  contained  in  Section  3.3.2. 
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3.3  Converting  array  time  series  data  to  traveltimes 


8.3.1  Identifying  the  arrivals 

Examples  of  the  data  that  were  received  at  the  hydrophone  array  were 
given  in  Figures  2.5a,  2  5b  and  2.5c.  Each  of  the  arrivals  represent  acoustic 
energy  that  has  traveled  along  a  different  path.  Figure  2.4  illustrates  a  few  of 
the  paths  that  sound  may  take  in  traveling  from  the  shot  to  the  array.  In  order 
to  utilize  these  data,  it  is  necessary  to  identify  the  raypath  associated  with  each 
arrival  on  the  array  time  series. 

To  associate  an  arrival  with  a  raypath,  time  domain  beamforming  is 
used  to  estimate  the  angle  of  arrival  for  each  arrival  at  the  array.  Once  this  is 
known,  a  ray  leaving  the  array  at  its  angle  of  arrival  may  be  traced.  A  velocity 
profile  for  the  water  column  was  derived  from  XBT  data  (to  1800  meters)  taken 
during  the  experiment.  Archival  data  was  used  to  extend  this  profile  from  1800 
meters  to  the  sea  floor.  By  appending  a  model  sedimentary  velocity  profile  [241 
to  the  velocity  profile  of  the  water  column,  bottom  interacting  rays  may  be 
traced  as  both  bottom  reflections  and  bottom  refractions.  This  is  illustrated  in 
Figure  3.1a.  Here  it  is  assumed  that  the  intersection  of  the  direct  path  (D)  and 
surface  reflection  (S)  gives  the  approximate  location  of  the  shot.  The  ray  that  is 
traced  as  a  bottom  reflection  (B)  misses  the  shot  point  by  a  wide  margin  while 
the  ray  traced  as  a  bottom  refraction  -  surface  reflection  (BR-S)  comes  very 
close.  The  conclusion  is  that  this  particular  arrival  refracted  in  the  sediment 
and  was  reflected  by  the  surface  before  reaching  the  array.  Figure  3.1b  shows 
all  of  the  raypaths  associated  with  arrivals  at  the  array  for  this  shot  It  should 
be  mentioned  at  this  point  that  the  true  velocity  profile  in  the  sediment  is  unk¬ 
nown  -  it  is  what  we  are  after  If  all  of  the  raytraces  of  bottom  interacting 
arrivals  had  come  as  close  as  shown  in  Figure  3.1b,  then  the  model  velocity 
profile  used  for  the  raytracing  would  be  an  excellent  estimate  and  there  would 
have  been  no  need  to  continue.  This  shot  was  particularly  good  and  was  chosen 
to  illustrate  the  procedure  for  identifying  the  arrivals.  Although  most  of  the 
other  shots  were  less  ideal,  there  was  enough  of  a  difference  in  the  raypaths  to 
clearly  identify  the  paths  they  had  taken  Angles  for  all  arrivals  are  given  in 
Table  3.1.  Angles  are  measured  from  the  horizontal;  arrivals  coming  from  below 
the  array  have  positive  angles,  those  coming  from  above  have  negative  angles. 


8. 8.2  Shot  locations 

A  constant  gradient  raytrace  program  was  used  to  estimate  shot  loca¬ 
tions.  Rays  were  traced  from  the  array  at  the  angle  of  arrival  for  the  direct 
water  path  (D)  until  range  along  the  ray  matched  the  radar  range  recorded 
when  the  charge  was  launched  The  depth  at  this  range  is  assumed  to  be  the 
shot  depth.  Times  to  that  range  and  depth  were  also  obtained  from  the  ray- 
trace. 
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(a) 


(b) 


Figure  3.1  Associating  an  arrival  at  the  array  with  a  raypath. 

(a)  A  ray  leaving  the  array  at  its  angle  of  arrival  is  traced  as  both  a  bottom  reflection  and  as  a 
bottom  refraction.  The  reflection  misses  the  shot  point  (assumed  to  be  the  intersection  of  the 
direct  and  surface  reflected  paths)  while  the  refraction  comes  very  close.  The  conclusion  is  that 
this  particular  arrival  refracted  in  the  sediment  and  was  reflected  by  the  surface  before  reaching 
the  array,  (b)  All  raypaths  for  a  medium-range  shot. 
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Anele  ol 

Arrival 

- - - - 

Shot 

Radar  Range 

(meters) 

D 

BR 

S 

BR-S 

S-BR 

S-BR-S 

38 

3704 

4.089 

34.166* 

-42.647 

-56  923* 

58.844* 

13 

3778 

1.022 

34.786* 

-40.822 

-57.871* 

58.517* 

1  _ 

39 

5741 

3.407 

29.751 

-30.934 

-44.527* 

_ 

15 

5834 

3  919 

29  164 

-31.532 

-44.289* 

_ 

40 

7408 

2.900 

21.620 

-24.769 

-37.958* 

_ 

17 

7871 

2.900 

19.981 

-23.278 

-36,250* 

37.742* 

-51.670* 

41 

9075 

2.896 

17.286 

-20.526 

-42.416 

_ 

-45.979* 

19 

9816 

3.066 

16.219 

-19.258 

-37.527 

35.410 

-43.580* 

37 

11575 

3.066 

13.949 

-16.601 

-28.387 

31.731 

-40.822* 

21 

11631 

3.237 

13.576 

-16.397 

-28.194 

29.359 

_ 

43 

12594 

3.407 

13.052 

-15  511 

-27.424 

29  751 

-44.527 

23 

13890 

3  237 

12.179 

-13.401 

-24.582 

26.089 

-39.706 

36 

14075 

3.407 

11  483 

-13.576 

-23.835 

25  333 

-38.391 

44 

14168 

3.578 

11  657 

-13  752 

-24  022 

25  711 

-38.391 

45 

15742 

3  919 

10  963 

-12  353 

-21.987 

23.278 

-31  133 

25 

16112 

4.089 

10  789 

-11  831 

-20.890 

22.539 

-32.536 

46 

16946 

4.260 

10.443 

-11.483 

-20.526 

21.255 

-31.332 

27 

18057 

4.602 

** 

-10.270 

-18.538 

19.439 

-28.775 

47 

18520 

4  943 

** 

-10.616 

-18.538 

19.078 

-29.359 

29 

20187 

5  456 

** 

-9  233 

-16.574 

17.822 

-26.089 

48 

20650 

5  542 

** 

-9.406 

-16.752 

17.464 

-26  469 

31 

22039 

5.798 

** 

-8  544 

-15.511 

16.574 

-24.208 

49 

22594 

6.141 

** 

-8.716 

-15.688 

16.042 

-24.208 

33 

23706 

6  483 

** 

-7.684 

-14.806 

15.865 

-22.908 

50 

24817 

7  000 

** 

-8.028 

-14.454 

14.982 

-22.723 

51 

26576 

7.513 

** 

-7.684 

-14.102 

13.567 

-21.072 

52 

28613 

7.856 

** 

-7  341 

-13.401 

12.179 

-19.439 

53 

30002 

8.028 

** 

-6  655 

-12.702 

11.310 

-18.001 

54 

32410 

8.028 

** 

-6  141 

-12.528 

10.616 

-16.930 

55 

34355 

- 

** 

-5.970 

-12.005 

10.443 

-16.042 

56 

36299 

- 

** 

-5.627 

-12.353 

9.751 

-15.865 

*  indicates  the  bottom  interaction  was  a  reflection  (not  a  refraction) 
**  indicates  path  does  not  exist 


Table  3.1  Angle  of  arrival  at  the  hydrophone  array.  Angles  are  measured 
from  the  horizontal;  arrivals  coming  from  below  the  array  have  positive  angles, 
those  coming  from  above  have  negative  angles. 
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Initially,  a  water  column  velocity  profile  derived  from  1800  meter  XBT 
data  taken  during  the  experiment,  and  archival  data  below  1800  meters  was 
used  to  trace  the  rays.  It  soon  became  evident  that  this  profile  was  inadequate. 
Many  of  the  rays  ended  up  far  from  the  assumed  shot  depth  of  1820  meters.  At 
short  ranges,  rays  refracted  too  much.  As  shot  range  increased,  rays  tended  to 
refract  too  little,  giving  greater  shot  depths  with  increasing  range.  The  steps 
involved  in  finding  velocity  and  attenuation  profiles  in  the  sediment  critically 
depend  upon  having  an  accurate  velocity  profile  in  the  water  column.  When 
using  the  direct  water  arrival  to  locate  the  shot,  the  ray  is  nearly  horizontal  and 
is  very  sensitive  to  gradients  in  the  water  column.  A  slight  increase  or  decrease 
in  bending  of  the  ray  translates  to  a  large  decrease  or  increase  in  the  final 
range.  An  accurate  velocity  profile  in  the  water  column  is  also  needed  when 
finding  sediment  traveltimes.  To  determine  the  range  and  time  that  a  ray 
spends  in  the  sediment,  the  range  and  time  that  the  ray  spends  in  the  water 
column  must  be  subtracted  out.  Sediment  refracted  arrivals  that  are  received 
at  the  array  spend  a  much  larger  amount  of  their  time  in  the  water  column 
than  in  the  sediment.  The  water  path  contribution  is  such  a  large  portion  of  the 
total  traveltime,  that  even  small  errors  in  removing  the  water  path  may  be  of 
the  same  order  as  the  traveltime  in  the  sediment. 

In  order  to  come  up  with  a  better  estimate,  the  direct  water  path  data 
were  inverted  to  find  a  velocity  profile  in  the  lower  portion  of  the  water  column 
(from  the  array  at  1500  meters,  to  the  sea  floor).  As  the  shot  range  increases, 
the  direct  path  rays  refract  deeper  as  shown  in  Figure  3.2.  The  data  set  pro^ 
vided  by  these  arrivals  is  of  the  same  form  as  that  which  is  used  to  find  the  velo¬ 
city  profile  in  the  sediment  (i.e.  T(X)  or  one  of  its  parameterizations,  X(p),  T(p), 
or  Ap))-  Inversion  of  the  data  is  a  straight-forward  problem  since  the  velocity 
profile  in  this  part  of  the  water  column  increases  monotonically  to  the  floor. 
The  least  squares  method  used  to  invert  the  data  is  described  in  detail  in  Sec¬ 
tions  3.4  and  3.6.1,  and  in  Appendix  Al. 

XBT  data  must  be  used  to  obtain  a  velocity  profile  in  the  upper  portion 
of  the  water  column  (from  the  surface  to  the  array).  When  looking  up  from  the 
array  toward  the  surface,  the  sound  channel  provides  a  low  velocity  zone  which 
makes  an  inversion  of  surface  reflected  data  a  difficult  problem.  Since  there  is 
much  scatter  in  the  XBT  data,  the  accuracy  of  the  velocity  profile  in  the  upper 
water  column  is  likely  to  be  poor.  Also,  since  the  experiment  was  carried  out 
over  a  two  day  period  with  shot  runs  being  conducted  between  noon  and  mid¬ 
night,  the  velocity  profile  in  the  upper  part  of  the  water  column  varied  as  a 

function  of  time.  At  best,  only  an  average  velocity  profile  for  this  region  could 
be  obtained. 

Although  a  few  shots  still  missed  the  nominal  depth  of  1820  meters  by 
a  large  margin,  overall  the  results  were  improved  using  the  velocity  profile 
obtained  from  the  above  method.  The  sound  speed  profile  in  the  water  column 
used  to  reduce  the  data  is  listed  in  Table  3.2  and  illustrated  in  Figure  3.3.  Shot 
locations  are  given  in  Table  3.3. 
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Figure  3.2  Raypaths  in  the  water  column 

Direct  paths  sample  deeper  parts  of  the  water  column  as  shot  range  increases. 

The  data  set  provided  by  these  arrivals  can  be  used  to  obtain  the  velocity  profile 
in  the  lower  portion  of  the  water  column. 


1480  1485  1490  1495  1500  1505 

velocity  (meters/second) 


Figure  3.3  Velocity  profile  in  the  water  column. 

Lower  portion  of  the  velocity  profile  was  obtained  from  inversion  of  direct 
water  path  data  (illustrated  in  Figure  3.2).  Upper  portion  of  the  profile  was 
obtained  from  XBT  data. 
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Velocity  Profile 

- 

Water  Column 

depth 

velocity 

depth 

velocity 

(meters) 

(meters/sec) 

(meters) 

(meters/sec) 

0. 

1494  0 

712.5 

1480.9 

16  1 

1492.4 

819.8 

1481.0 

37.5 

1490.8 

870  8 

1481.1 

53.9 

1489.9 

900.0 

1481.2 

75.0 

1489.0 

937.5 

1481.3 

109.2 

1487.8 

975.0 

1481.4 

150.0 

1486.6 

1016.3 

1481.5 

185  0 

1485.7 

1050.0 

1481.6 

225.0 

1484.8 

1136.3 

1482.0 

259.8 

1484.1 

12000 

1482  4 

300.0 

1483  4 

1358.9 

1483.8 

333.7 

1482  9 

1500.0 

1485.2 

375.0 

1482  4 

1555.0 

1485.7 

409.8 

1482  1 

1831.0 

1488.0 

450  0 

1481.8 

2013.0 

1490.2 

483.6 

1481.5 

2172.0 

1492.5 

525.0 

1481.2 

2461.7 

1496.9 

600.0 

1481.0 

2800.0 

1502.0 

628.6 

675.0 

1480  9 
1480.8 

2853.0 

1502.8 

Table  3.2  Velocity  profile  in  the  water  column 


8.S.S  Sea  Floor  Depth 

Two  sources  of  information  for  estimating  the  sea  floor  depth  were 
used.  An  iterative  method  that  uses  data  collected  from  the  streamer  hydro¬ 
phone  was  used  to  estimate  both  shot  and  sea  floor  depths.  Not  all  shots  had 
streamer  data  that  could  be  used,  so  the  number  of  estimates  using  this  method 
is  limited. 

A  second  estimate  of  sea  floor  depth  was  obtained  from  an  oceano¬ 
graphic  map  of  the  Monterey  fan  (‘Oceanographic  Data  of  the  Monterey  Deep 
Sea  Fan’,  1st  Edition,  June  1975  by  T.E.  Chase,  W.R.  Normark,  and  P.  Wilde. 
Prepared  at  the  Geologic  Data  Center,  Scripps  Institution  of  Oceanography, 
IMR  Tech  Rept.  Series  TR-58).  Both  sources  indicate  a  slope  of  about  .4  degrees 
out  to  a  range  of  20  km,  and  a  slope  of  1.4  degrees  beyond  that  range.  Sea  floor 
depth  as  a  function  of  range  from  FLIP  is  given  below. 
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Shot  Location 


Shot 


Radar 

Range 


Arrival 
Angle  (D) 


Raytrace  along  direct  path 


Range 


Depth 


Time 


39 

15 

40 
17 

41 
19 
37 
21 

43 
23 
36 

44 

45 
25 

46 
27 

47 
29 

48 
31 

49 
33 

50 

51 

52 


5741 

5834 

7408 

7871 

9075 

9816 

11575 

11631 

12594 

13890 

14075 

14168 

15742 

16112 

16946 

18057 

18520 

20187 

20650 

22039 

22594 

23706 

24817 

26576 

28613 


3.407 
3  919 
2.900 

2  900 
2.896 
3.066 
3.066 
3.237 
3.407 
3.237 

3  407 
3  578 

3  919 
4.089 
4.260 

4  602 

4  943 

5  456 
5.542 

5  798 

6  141 
6.483 

7  000 
7.513 
7  856 


5744 

5826 

7424 

7880 

9080 

9868 

11534 

11657 

12581 

13876 

14075 

14164 

15741 

16119 

16955 

18054 

18533 

20182 

20649 

22038 

22592 

23709 

24817 

26572 

28610 


1799 

1853 

1771 

1775 
1778 
1805 

1794 
1820 
1817 
1778 

1776 

1795 
1794 
1813 
1812 
1797 
1820 
1819 
1794 
1729 
1768 
1760 

1817 

1818 
1706 


3.867 

3.923 

4.996 

5.303 

6.109 

6.639 

7.759 

7.841 

8  462 
9.333 

9  467 
9.526 

10.585 

10.839 

11.400 

12  138 
12.458 

13  564 
13  878 
14.812 
15.180 
15.928 
16.665 
17.836 
19.205 


Table  3.3  Shot  locations  Rays  were  traced  from  the  array  at  the  angle  of  ar¬ 
rival  along  the  direct  water  path  (D)  until  range  along  the  ray  matched  the  ra¬ 
dar  range  recorded  when  the  shot  was  launched.  The  depth  at  this  range  is  as¬ 
sumed  to  be  the  shot  depth.  Ranges  and  depths  are  reported  in  meters. 


8.8.4  Removing  the  water  path  contribution 

Having  identified  the  arrivals  and  found  the  shot  locations,  rays  may  be 
traced  from  the  shot  and  from  the  array  down  to  the  sea  floor.  A  time  origin 
was  not  available  so  times  in  the  sediment  were  found  from  the  time  difference 
between  the  direct  and  refracted  paths. 
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Sea  Floor  Depth 

range  from 

depth 

range  from 

depth 

FLIP  (meters) 

(meters) 

FLIP  (meters) 

(meters) 

0 

2800 

9700 

2861 

2100 

2813 

11400 

2871 

3100 

2819 

11630 

2873 

3900 

2824 

20190 

2932 

4900 

2831 

22040 

2973 

5300 

2833 

24820 

3009 

6900 

2843 

26580 

3086 

8000 

2850 

32410 

3232 

9100 

2857 

Table  3.4  Sea  floor  depth 


water 

column 


Figure  3.4  Removing  the  water  path  contribution. 
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With  reference  to  Figure  3.4,  define: 


Tr  Total  time  for  sediment  refracted  raypath  (including  time  spent 
in  the  water  column) 

TD  Total  time  for  the  direct  path  (water  column) 

<1  Time  from  the  array  to  the  floor 
<2  Time  from  the  shot  to  the  floor 
ts  Time  in  the  sediment 

The  total  time  for  the  sediment  refracted  raypath  is 
Tr  =  +  <2  +  h 

The  time  difference  between  the  sediment  and  direct  path  is 
At  =  Tr  -Td 

=  +  £2  +  ~  Tp 

Therefore  the  time  that  the  ray  spends  in  the  sediment  is  given  by 
tR  —  At  Tp  —  —  <2 

All  of  the  variables  on  the  right  side  of  (3.1)  may  be  found. 

At  is  found  from  the  array  time  series 

(time  difference  between  the  direct  and  sediment  refracted  paths). 

Tp  is  obtained  by  raytracing  (or  streamer  data  if  available). 

t\t2  are  obtained  by  raytracing. 


Ranges  in  the  sediment  corresponding  to  times  ts  are  found  by  raytrac¬ 
ing.  The  ranges  associated  with  times  and  t2  are  subtracted  from  the  shot 
range  (from  FLIP  to  the  shot).  Angles  used  to  trace  rays  from  the  array  are 
listed  in  Table  3.1.  Shot  ranges  are  given  in  Table  3.3.  Time  differences 
between  the  direct  and  refracted  paths  (At)  are  listed  in  Table  3.5. 

Data  from  longer  range  shots  contain  a  considerable  amount  of  error. 
So  much  of  the  ray’s  total  range  and  time  is  spent  in  the  water  column  that 
errors  can  be  larger  than  the  range  and  time  spent  in  the  sediment.  For  this 
reason,  shots  beyond  a  range  of  20  km  are  not  used. 
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Arrival  Time  Differences 

difference  in  time  fms)  between  direct  and  refracted  paths 

Shot 

BR 

BR-S 

S-BR 

S-BR-S 

39 

335 

15 

321 

40 

242 

17 

221 

41 

184 

1246 

19 

158 

1102 

1260 

37 

127 

922 

1054 

21 

123 

891 

1042 

43 

109 

832 

974 

2537 

23 

98 

769 

845 

2007 

36 

90 

740 

831 

1992 

44 

89 

745 

845 

1940 

45 

76 

680 

757 

1793 

25 

63 

645 

745 

1744 

46 

61 

627 

708 

1675 

27 

50 

582 

652 

1532 

47 

569 

656 

1530 

29 

520 

600 

1404 

48 

519 

591 

1389 

31 

479 

549 

1291 

49 

467 

547 

1281 

33 

450 

515 

1211 

50 

439 

508 

1184 

51 

421 

477 

1117 

52 

405 

464 

1056 

Table  3.5  Arrival  time  differences  between  direct  and  refracted  paths. 
Differences  in  arrival  times  were  found  from  hydrophone  15  of  the  array  time 
series.  This  is  roughly  the  middle  of  the  upper  half  (high  gain  portion)  of  the  ar¬ 
ray.  Arrival  angles  given  in  Table  3.1  were  found  for  this  section  of  the  array. 


3.4  Re-parameterization  of  the  data 

T(X)  is  not  a  particularly  good  form  for  the  data  when  trying  to  infer 
velocity  structure.  Small  changes  in  the  traveltime  curve  cause  large  changes  in 
the  model.  It  is  common  to  parameterize  range  and  time  as  functions  of  ray 
parameter  p.  While  this  results  in  a  more  stable  form  of  the  data,  X(p)  and 

T(p)  are  sensitive  to  errors  made  in  estimating  p,  and  to  lateral  inhomo¬ 
geneities. 

A  linear  combination  of  X(p)  and  T{P),  called  the  delay  time,  is  defined 
as.  r(p)  =  T(p)  -pX(p).  Unlike  X(p)  and  T(p),  the  delay  time  is  a  monotonically 
decreasing  function  and  is  sensitive  only  to  second  order,  to  errors  in  p  [10], 


31 


R.  K.  Brienzo 


Suppose  there  is  an  ( X,T )  pair  which  corresponds  to  a  true  ray  parameter  of  p0. 
If  an  incorrect  ray  parameter  px  is  estimated  for  this  data  pair,  errors  will  be 
introduced  into  the  X(p)  and  T(p)  data  since  the  point  (X(Pi),  T(Pl ))  is  off  the 
correct  traveltime  curve.  t(p )  is  somewhat  insensitive  to  these  errors;  even 
though  an  incorrect  value  of  p  =p,  is  used,  t(Pi)  will  still  be  correct  to  first 
order. 

The  price  that  is  paid  for  this  insensitivity  to  errors  is  a  loss  of  some 
detail  in  the  velocity  profile.  The  smoothness  of  the  r(p)  curve  can  make  it 
difficult  to  recover  details  in  the  velocity  profile.  Orcutt,  MacKenzie  and 
McClain  clearly  illustrate  this  [53,  Fig.  2],  r(p)  and  X(p)  data  are  generated  for 
two  distinctly  different  velocity  profiles.  While  the  two  X(p)  curves  are  quite 
different,  the  two  t(p)  curves  are  very  similar.  Given  noise  in  the  data,  it  would 
be  difficult  to  recover  the  two  velocity  profiles  from  the  j(p)  data.  Although  the 
two  velocity  profiles  could  be  recovered  from  the  X(p)  data,  this  data  is  sensitive 
to  errors  in  estimating  p  and  to  lateral  inhomogeneities. 

Expressing  X  as  a  function  of  p  is  often  a  difficult  step.  By  definition, 
the  ray  parameter  is  given  by  p  =  {dT/dX),  the  slope  of  the  traveltime  curve  at 
a  particular  (X,T)  pair.  Phis  seldom  produces  acceptable  results  however,  since 
numerical  differentiation  of  noisy  data  is  an  unstable  procedure. 

A  standard  technique  for  estimating  p  is  to  fit  a  polynomial  through 
T{X)  and  then  differentiate  the  polynomial.  This  is  better  than  trying  to  find 
the  slope  from  discrete  points  along  the  traveltime  curve,  but  fitting  a  polyno¬ 
mial  to  the  data  can  smooth  over  small  variations  in  the  traveltime  curve  which 
provide  detail  in  the  velocity  profile.  Although  methods  exist  for  obtaining  suit¬ 
able  bounds  on  r(p)  without  actually  estimating  p  [34,  35,  52,  53],  converting  the 
data  into  a  function  of  p  is  a  difficult  and  error-prone  step  in  most  seismic  stu- 
dies. 

The  problem  of  having  to  estimate  p  from  the  traveltime  data  is  essen¬ 
tially  side-stepped  in  this  study.  Since  an  array  was  used  to  receive  the  acoustic 
energy,  the  ray  parameter  can  be  calculated  from  Snell’s  law.  If  a  medium  is 
laterally  homogeneous  and  horizontally  layered,  then  p  is  constant  along  the 
entire  raypath  and  p  =u  sin(9).  u(z),  the  reciprocal  of  velocity,  is  called  slow¬ 
ness:  u{z)  =1  /v(z).  e  is  the  angle  the  ray  forms  with  the  vertical.  The  angle  of 
arrival  may  be  found  by  beamforming  and  thus  p  can  be  calculated.  (Note: 
angle  of  arrival  values  reported  in  Table  3.1  were  measured  with  respect  to  the 
the  horizontal,  so  6  =  90°- angle  of  arrival.) 

X(p)  will  be  used  in  the  least  squares  inversion  in  order  to  gain  the 
detail  available  from  this  form.  Having  an  good  estimate  of  p  associated  with 
each  X  eliminates  the  necessity  of  fitting  a  polynomial  through  the  traveltime 
data,  or  using  r(p).  Details  that  might  otherwise  be  lost  will  thus  be  preserved. 
r(p)  will  be  used  in  the  linear  programming  inversion  to  obtain  extremal  bounds. 
Here  the  intent  is  to  produce  bounds,  the  smoothing  that  occurs  from  using  r(p) 
is  perfectly  acceptable.  X(p)  data  will  be  incorporated  into  this  method  as  addi¬ 
tional  constraints. 

The  traveltime  data  and  its  various  parameterizations  are  listed  in 
Table  3.6  and  illustrated  in  Figures  3.5a-d. 
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Sediment  traveltime  data 

p  X  103 

X 

T 

r{p) 

.49692 

3116 

1.933 

.38460 

.51783 

3033 

1.860 

.28942 

.52756 

2713 

1.708 

.27673 

.52756 

2778 

1.694 

.22844 

.53380 

2896 

1.807 

.26112 

.56745 

1918 

1.253 

.16463 

.57248 

1997 

1.285 

.14176 

.57493 

2059 

1.345 

.16122 

.58436 

1665 

1.074 

.10104 

.58436 

2209 

1.417 

.12615 

.58664 

2387 

1.529 

.12869 

.58776 

1742 

1.131 

.10712 

.58997 

1549 

1.006 

.09214 

.59215 

1366 

.927 

.11812 

.59322 

1447 

.958 

.09961 

.59744 

2019 

1.309 

.10277 

.60450 

1590 

1.047 

.08585 

.60450 

1820 

1.168 

.06781 

.60644 

1865 

1.208 

.07699 

.60836 

1604 

1.037 

.06119 

.61208 

1753 

1.157 

.08402 

.61479 

1739 

1.145 

.07588 

.61568 

1496 

.979 

.05794 

.61829 

1945 

1.261 

.05843 

.62167 

1751 

1.143 

.05446 

.62413 

1967 

1.300 

.07234 

.62573 

1320 

.867 

.04104 

.62730 

1606 

1.054 

.04656 

.62884 

1547 

1.021 

.04819 

.63035 

2072 

1.367 

.06092 

.63257 

1176 

.775 

.03110 

.63472 

1111 

.733 

.02783 

.63611 

1166 

.782 

.04030 

.63816 

1264 

.842 

.03537 

.63816 

1813 

1.199 

.04202 

.64078 

1426 

.943 

.02925 

.64268 

1092 

.725 

.02319 

.64512 

1140 

.762 

.02656 

.64630 

1289 

.852 

.01892 

.65323 

1048 

.697 

.01242 

.65428 

880 

.588 

.01223 

.65569 

1205 

.802 

.01189 

65793 

1078 

.722 

.01275 

.65920 

572 

.384 

.00694 

.66080 

704 

.476 

.01080 

.66119 

745 

.495 

.00241 

.66193 

598 

.403 

.00717 

Table  3.6  X,  T,  and  r  data 
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Travel  time  data  T(X) 


Figure  3.5  (a)  T(X)data,  (b)  T(p)data. 
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3.5  Sources  of  error 

There  are  several  sources  of  error  responsible  for  scatter  in  the  travel¬ 
time  data.  These  are  present  from  uncertainties  occurring  during  the  experi¬ 
ment  and  also  errors  introduced  during  processing  of  the  data.  Each  source  of 
error  will  be  discussed  and,  when  possible,  an  estimate  made  of  its  effect  on  the 
data.  In  most  cases,  only  an  order  of  magnitude  estimate  can  be  made. 


The  water  column  and  sediment  are  assumed  to  consist  of  horizontal, 
laterally  homogeneous  layers.  Velocity  is  assumed  to  increase  with  depth.  Vari¬ 
ations  in  topography  (such  as  sloping  interfaces  or  lateral  inhomogeneities)  and 
velocity  anisotropy  violate  the  model  assumptions  and  introduce  errors.  The 
region  in  which  the  experiment  was  conducted  is  relatively  flat,  but  in  an  area  of 
rapid  deposition  and  turbidity  currents,  lateral  inhomogeneities  and  velocity 
anisotropv  are  to  be  expected.  Alternating  layers  of  high  and  low  velocity 
material  where  lithification  begins  to  occur  can  also  cause  scatter  in  the  travel- 
time  data.  Without  more  information  about  subsurface  layering  it  is  difficult  to 
assess  the  magnitude  of  the  errors  introduced  by  these  effects. 


Arrival  time  errors. 

Arrival  time  errors  occur  when  arrival  times  are  picked  from  the  array 
time  series.  These  errors  are  relatively  small,  typically  on  the  order  of  a  few 
milliseconds  Even  for  the  deeply  refracted  rays  whose  onset  time  is  less  dis¬ 
tinct,  uncertainties  in  arrival  time  picks  are  at  most  about  20  milliseconds.  For 
a  sound  speed  of  1500  meters/sec  this  corresponds  to  an  error  in  range  of  at 
most  30  meters  (usually  much  less). 

Errors  m  X  and  T  that  show  up  from  errors  in  estimating  p. 

Dorman  10  has  shown  that  errors  in  estimating  p  can  have  a 
significant  effect  on  uncertainties  in  X.  In  parameterizing  the  data  as  X(p),  ray 
parameter  p  is  the  independent  variable  and  range  X  is  the  dependent  variable. 
The  independent  variable  is  usually  considered  to  be  free  of  error,"  all  errors 
show  up  in  the  dependent  variable. 

The  total  error  in  X  is  made  up  of  two  terms. 

Et  =  Ex  +  Ep 

There  are  measurement  errors  in  X  (denoted  Ex).  These  would  arise  from  errors 
affecting  range  estimates,  such  as  incorrect  shot  depth  or  radar  ranges,  a  less 
than  perfect  water  column  velocity  profile,  lateral  inhomogeneities,  etc.’  There 
are  also  errors  in  X  that  show  up  due  to  incorrect  estimates  of  the  ray  parame¬ 
ter  p  (denoted  Ep). 

Assume  there  is  a  traveltime  pair  (X,  T)  that  corresponds  to  ray  param¬ 
eter  p  =p0,  but  an  incorrect  estimate  of  p  =Pl  is  made.  The  estimate  X(Pl)  does 
not  lie  on  the  true  traveltime  curve  (the  correct  value  is  X(Po)  ).  The  error  in  X 
due  to  the  incorrect  value  of  p  is 

Ep  =  X(Pl)  —  X(Po)  (3.2) 
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Let  the  error  in  estimating  p  be  defined  as  Sp  =  Pl  -Po,  and  express  X(Pl)  by  a 
Taylor  series  expansion  about  X(Po). 

X(Pi)  =  X(Po)  +  i'^')poSP  +  0(Sp2)  + 


Dropping  the  higher  order  terms  and  substituting  into  (3.2)  gives 


(3.3) 


Since  p 


dT 
dX  ; 


_dp_ 

dX 


d2T 

dX2 


and  we  have  that 


dX  _  (  dp  W  _  /  d2Tyi 
dp  V  dX}  ~  \  dX2  ’ 


(3.4) 


Using  equation  (3.4)  in  (3.3)  gives 


Ev  =  ( 


£T_y 1 
dX2 ’ 


Sp 


(3.5) 


d  2  y 

Note  that  ^he  curvature  of  the  traveltime  curve.  As  the  traveltime  curve 

approaches  a  straight  line,  the  reciprocal  of  the  curvature  becomes  infinite  fa 
straight  line  has  =-0  ).  Equation  3.5  shows  that  errors  in  estimating  p  can 

be  greatly  magnified.  It  is  possible  that  errors  due  to  uncertainties  in  estimating 
P  ( Ep )  could  dominate  errors  due  to  uncertainties  in  range  measurements  (Ex). 

To  examine  the  effect  of  uncertainties  in  p  in  the  context  of  this  data 
set,  the  curvature  of  the  traveltime  data  will  be  estimated  and  used  with  a  rea¬ 
sonable  value  of  Sp  to  find  the  resulting  error  in  range.  Fitting  a  quadratic 
through  the  traveltime  data  (Table  3.6)  yields  the  equation 

T  -  .0012  +  T^XIO^X  -  8.5  XlO-8.^2 

so  that 

d2T\~t 

~Xp2J  w  ®  X  108  meters2 / sec 

From  equation  3.5  the  error  in  X  due  to  errors  in  estimating  p  is 

EP  =  ("fjr)  SP  =  (6  X  108)  Sp  (3.6) 


Since  the  data  were  received  at  an  array,  p  is  calculated  from  Snell’s 
law:  p  =  ustn(ff).  Given  the  high  signal  to  noise  ratio  of  the  data,  most  of  the 
error  in  P  should  be  due  to  the  finite  number  of  look  directions  that  were  avail¬ 
able  (the  beamformer  used  for  angle  of  arrival  estimation  only  allowed  an 
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integer  number  of  sample  delays).  This  error  in  p  for  the  steepest  ray  (where 
the  effect  is  most  severe)  amounts  to  2X10-8  sec /meter .  This  results  in  an  error 
Ep  on  the  order  of  10  meters.  Having  an  array  and  high  signal  to  noise  ratio 
data  has  kept  this  type  of  error  small. 


Errors  in  removing  the  water  path  contribution. 

The  largest  errors  occur  when  the  range  and  time  spent  in  the  water 
column  are  subtracted  from  the  ray’s  total  range  and  time.  The  time  spent  in 
the  water  column  is  much  greater  than  the  time  spent  in  the  sediment,  so  even 
small  errors  in  subtracting  the  water  path  are  a  large  percentage  of  the  sedi¬ 
ment  traveltime.  These  errors  are  primarily  due  to  having  an  inexact  velocity 
profile  for  the  water  column  This  affects  both  the  actual  tracing  of  rays  to  the 
sea  floor  as  well  as  the  determination  of  shot  locations  (see  Section  3  3.2).  An 
additional  error  in  shot  location  occurs  if  the  explosive  charge  does  not  descend 
vertically  (it  may  drift  horizontally  due  to  currents).  Other  errors  include  those 
due  to  a  finite  number  of  look  directions  (about  25  meters)  and  uncertainties  in 
radar  ranges,  which  should  be  accurate  to  within  30  meters.  Since  the  correct 
velocity  profile  is  not  available,  it  is  difficult  to  know  the  error  introduced  by 
using  an  estimate. 


Total  error 

An  examination  of  scatter  in  the  data  indicates  a  total  uncertainty  in 
range  of  about  one  kilometer,  and  in  time  of  about  6  seconds.  Most  of  this 
occurs  in  removing  the  water  path.  The  biggest  reduction  in  scatter  of  the  data 
would  result  from  an  improved  estimate  of  the  velocity  profile  in  the  water 
column. 


3.6  Inverting  the  data 


Two  methods  will  be  used  to  invert  the  data.  The  problem  will  first  be 
solved  using  a  least  squares  approach.  This  solution  yields  velocity  as  a  function 
of  depth.  A  second  inversion,  using  a  linear  programming  method,  produces 
bounds  on  the  velocity  profile.  All  models  that  are  consistent  with  the  data  lie 
within  these  bounds.  These  bounds  indicate  in  an  unambiguous  way,  what  can 
be  inferred  from  a  finite  set  of  data  with  errors.  A  single  velocity  profile  is 
needed  for  acoustic  modeling;  the  least  squares  solution  provides  this  The 
linear  programming  solution  provides  a  way  to  assess  the  quality  of  the  solution. 
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8.6.1  Least  squares  solution 

In  this  approach  the  sediment  is  modeled  as  a  horizontally  layered, 
laterally  homogeneous  medium.  Velocity  is  assumed  to  increase  with  depth! 
Layers  in  the  model  are  defined  by  specifying  a  slowness  value  for  each  interface 
(depths  to  the  interfaces  are  not  known).  Each  layer  is  described  by  a  constant 
slowness  gradient.  The  least  squares  solution  yields  the  unknown  slowness  gra¬ 
dients.  Layer  thicknesses  may  be  found  from  the  gradients  and  the  interface 
slowness  values  ( dz/du  X  du  =  dz ).  Velocity  as  a  function  of  depth  is  then 
found  by  integrating  the  layer  thicknesses.  The  matrix  equation  that  is  inverted 
to  obtain  the  unknown  gradients  is  derived  in  Appendix  Al.  Details  of  the 
model  are  also  described  in  this  appendix. 

A  singular  value  decomposition  15,  39,  41,  66]  was  used  to  invert  the 
matrix  equation  (A1.12).  No  zero  singular  values  were  obtained  for  any  of  the 
models  that  were  examined.  The  largest  difference  in  eigenvalues,  occurring  for 
a  five  layer  model,  was  Xmax/Xmin  =  6.34,  indicating  a  well  conditioned  matrix.  In 
all  cases  that  will  be  presented,  solutions  were  constructed  using  the  maximum 
number  of  eigenvectors.  The  maximum  number  of  eigenvectors  is  equal  to  the 
number  of  layers  in  the  model  (the  number  of  columns). 

The  velocity  profile  in  the  upper  sediments  is  usually  represented  by  a 
monotonically  decreasing  function  that  changes  smoothly  with  depth  [22,  23,  24 
Defining  constant  gradients  in  each  layer  of  the  model  allows  representation  of 
such  a  function  by  only  a  few  layers.  At  best,  it  is  unnecessary  to  use  many 
more  layers,  and  at  worst,  is  likely  to  cause  problems  in  the  least  squares  inver¬ 
sion.  If  the  number  of  layers  in  the  model  approaches  or  exceeds  the  number  of 
data,  then  gradients  are  found  for  layers  in  which  few,  if  any  rays  turn.  This 
often  leads  to  physically  unreasonable  solutions. 

To  examine  the  effect  of  model  selection  on  the  solution,  X(p)  data  are 
inverted  for  models  ranging  from  one  to  five  layers  (Figures  3.6-3  10).  Three 
plots  for  each  inversion  will  be  presented.  The  first  plot  (a)  shows  the  velocity 
profile  obtained  from  the  inversion.  The  profile  suggested  by  Hamilton  is 
included  for  comparison.  In  (6)  and  (c),  theoretical  X(p)  and  t(p)  curves  are  gen¬ 
erated  from  the  velocity  profile  and  plotted  along  with  the  experimental  data. 

Figure  3.6  shows  the  result  for  a  one  layer  model.  Velocity  gradients  in 
the  upper  sediments  change  with  depth  (overall,  they  decrease  as  depth 
increases)  so  a  single  layer  would  not  be  expected  to  provide  a  good  fit.  As  seen 
in  Figure  3.6b,  the  the  X(p)  data  generated  from  the  velocity  profile  does  not  fit 
the  experimental  data  very  well  for  values  of  p  less  than  .55.  Most  of  the  data 
is  from  values  of  p  greater  than  .55,  so  the  emphasis  in  the  solution  is  focused 
more  on  this  region.  The  remaining  models  fit  the  data  more  closely.  As  sug¬ 
gested  in  Section  3.4,  there  is  not  much  difference  in  the  t (p)  curves  generated 
from  velocity  profiles  that  clearly  differ  Figure  3.11  shows  the  3  and  4  layer 
solutions  with  their  90%  confidence  bounds.  Details  of  these  velocity  profiles  are 
summarized  in  Tables  3.7  and  3.8. 

With  a  finite  data  set  there  are  an  infinite  number  of  models  that  will 
fit  the  data,  and  a  single  model  must  be  (often  somewhat  arbitrarily)  selected 
A  four  layer  model  adequately  describes  the  velocity  profile  without  running  into 
problems  with  data  density,  and  will  be  used  in  the  subsequent  analysis. 
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Figure  3.6  Inversion  using  a  one  layer  model. 

(a)  velocity  profile  obtained  from  inversion,  (b)  and  (c)  show  data  generated 
from  the  velocity  profile  in  part  (a)  (solid  line)  along  with  the  experimental  data. 
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Figure  3.7  Inversion  using  a  two  layer  model. 

(a)  velocity  profile  obtained  from  inversion,  (b)  and  (c)  show  data  generated 
from  the  velocity  profile  in  part  (a)  (solid  line)  along  with  the  experimental  data. 
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Figure  3.8  Inversion  using  a  three  layer  model. 

(a)  velocity  profile  obtained  from  inversion,  (b)and(c)  show  data  generated 
from  the  velocity  profile  in  part  (a)  (solid  line)  along  with  the  experimental  data. 
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Figure  3.10  Inversion  using  a  five  layer  model 

(a)  velocity  profile  obtained  from  inversion,  (b)  and  (c)  show  data  generated 
from  the  velocity  profile  in  part  (a)  (solid  line)  along  with  the  experimental  data. 
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Figure  3.11  90%  confidence  bounds  for  the  3  and  4  layer  models. 
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Velocity  Profil 

e  -  3  layer  model 

interface 

depth 

velocity 

gradient 

( meters ) 

(meters  /  sec} 

(sec-1) 

1 

0. 

1484.8 

2 

190  8 

1627.0 

.745 

3 

349  9 

1799.3 

1.083 

4 

620.7 

2012.4 

.787 

Table  3.7  Velocity  profile  in  the  sediment  for  a  three  layer  model. 


Velocity  Profile  -  4  layer  model 

interface 

depth 
( meters ) 

velocity 
(meters  /sec) 

gradient 

(sec-1) 

i 

0. 

1484.8 

2 

144  1 

1588  9 

.722 

3 

269.7 

1708.8 

.955 

4 

408.3 

1848.2 

1.006 

5 

621.6 

2012.4 

.770 

Table  3.8  Velocity  profile  in  the  sediment  for  a  four  layer  model. 


3.6.2  Linear  programming  solution 

The  solution  to  a  linear  programming  problem  minimizes  or  maximizes 
an  objective  function  subject  to  a  set  of  constraints.  In  this  case,  bounds  on  the 
input  data  (r(p),  X(p),  etc.)  form  the  set  of  constraints.  The  objective  function 
is  depth  to  a  particular  velocity.  Extremal  bounds  on  the  velocity  profile  are 
found  by  selecting  a  number  of  velocities  to  evaluate.  At  each  velocity  two 
linear  programming  problems  are  solved.  The  minimum  depth  at  that  velocity 
that  satisfies  the  constraints  is  found.  Next,  the  maximum  depth  at  that  velo¬ 
city,  subject  to  the  constraints,  is  found.  When  this  procedure  is  performed  at  a 
number  of  different  velocities,  extremal  bounds  are  swept  out.  These  bounds 
define  the  minimum  and  maximum  depths  possible  at  each  velocity  that  are  con¬ 
sistent  with  the  constraints  (bounds)  on  the  input  data  [14], 

One  advantage  of  this  method  is  that  additional  information  or  con¬ 
straints  may  be  easily  incorporated  into  the  problem.  Adding  constraints  on 
V(p)  data  restricts  the  slope  variations  possible  for  the  r(p)  data  and  can  tighten 
the  extremal  bounds  considerably. 

Input  bounds  on  the  r[p)  and  X(p)  data  are  shown  in  Figure  3.12. 
Extremal  bounds  found  using  only  r(p)  data  are  given  in  Figure  3.13a.  Much 
tighter  bounds  result  (Figure  3.13b)  when  the  X{p)  constraints  of  Figure  3.12b 
are  included. 
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Figure  3.12  Input  bounds  on  x  (p)  and  X(p)  data. 
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Figure  3.13  Linear  programming  solution. 

Extremal  bounds  (bold  lines)  are  plotted  along  with  the  least  squares  solution, 
(a)  inversion  of  t  (p)  data  only,  (b)  inversion  of  both  z  (p)  and  X(p)  data. 
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3.7  Forward  modeling 


To  provide  a  qualitative  check  on  the  results  of  Section  3,6.1,  theoreti¬ 
cal  curves  generated  from  the  velocity  profile  (Table  3.8)  will  be  compared  with 
the  experimental  data.  These  are  illustrated  in  Figures  3  14a-d. 

The  theoretical  curves  fit  the  data  reasonably  well  except  for  the  few 
data  at  the  largest  values  of  p.  Large  values  of  p  correspond  to  the  shallowest 
rays.  They  spend  nearly  all  of  their  time  in  the  water  and  merely  graze  the  sed¬ 
iment.  These  data  are  most  prone  to  error  since  such  a  large  portion  of  their 
traveltime  is  subtracted  off.  If  the  initial  velocity  used  in  the  inversion  is 
increased  from  1484.8  to  1494  m/sec,  the  resulting  velocity  profile  has  a  shal¬ 
lower  gradient  in  the  first  layer.  The  theoretical  curves  from  this  velocity  profile 
fit  the  large  p  data  better,  but  result  in  a  poorer  overall  fit  to  the  experimental 
data. 


3.8  Discussion 


Values  of  the  velocity  gradients  presented  in  Tables  3.7  and  3.8  are 
consistent  with  previously  published  values  for  this  type  of  sediment  [11,  22,  24, 
29],  The  way  these  gradients  change  with  depth  -  initially  shallow,  increasing  to 
about  one,  and  then  decreasing  -  differs  from  the  monotonically  decreasing  gra¬ 
dients  usually  reported.  This  feature  appeared  in  all  models  that  were  exam¬ 
ined,  and  velocity  profiles  which  do  not  exhibit  a  steeper  gradient  between  150 
and  400  meters  do  not  fit  the  experimental  data  very  well.  In  some  cases  this 
difference  may  be  due  to  the  way  in  which  the  data  were  analyzed.  Two  exam¬ 
ples  are  given  below. 

Velocity  profiles  obtained  from  reflection  data  are  found  by  averaging 
interval  velocities.  A  second  or  third  order  regression  through  these  points  is 
used  to  describe  velocity  as  a  function  of  depth.  A  smooth  profile  is  a  natural 
result  of  this  method. 

In  most  refraction  studies  it  has  been  necessary  to  fit  the  traveltime 
curve  with  a  polynomial,  and  differentiate  the  polynomial  to  convert  the  data 
into  a  function  of  p.  If  a  low  order  polynomial  is  fit  through  the  traveltime  data, 
small  details  which  give  structure  to  the  velocity  profile  may  be  lost.  Use  of  r(p) 
data  has  a  similar  smoothing  effect.  These  would  both  lead  to  smoother  velocity 
profiles.  In  this  study,  beamforming  at  the  array  allowed  the  direct  calculation 
of  p  and  permitted  the  use  of  V(p)  data.  Details  contained  in  the  traveltime 
data  were  preserved. 

An  increase  in  the  gradient  at  around  150  meters,  and  the  subsequent 
decrease  at  400  meters  indicates  that  a  change  in  some  physical  property  occurs 
near  those  depths.  There  are  several  physical  processes  that  might  account  for 
this  velocity  structure. 
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Figure  3.14  Forward  modeling. 

A  theoretical  curve  generated  from  the  velocity  profile  of  Table  3.8  is  plotted  along  with 
the  experimental  data,  (a)  T(X)data,  (b)  T(p)data. 
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Subsurface  topography  of  deep  sea  fans  can  be  very  complex.  Major 
events  (rise  or  fall  of  sea  level,  turbidite  flows,  etc.),  and  changes  in  channel 
direction  or  activity  can  create  irregular  structures  or  cause  primary  areas  of 
deposition  to  shift  with  time.  The  upper  portion  of  the  Monterey  fan  has  a  com¬ 
plicated  history  of  development  [27,  47,  48],  An  irregular  change  in  the  velocity 
gradients  with  depth  could  be  related  to  changes  in  sedimentation  history.  Velo¬ 
city  would  be  expected  to  change  as  age  and  sediment  type  vary. 

Mineralogical  changes  that  occur  with  depth  could  also  produce  a 
higher  gradient  at  mid-depth.  Lithification  of  the  unconsolidated  sediments  is 
thought  to  begin  between  200-400  meters  depth.  Cementation  of  the  mineral 
particles  increases  frame  rigidity.  This  higher  velocity  material  would  cause  a 
steeper  gradient.  After  lithification,  there  is  little  further  increase  in  frame  rigi- 
dity,  or  further  decrease  in  porosity.  Velocity  increases  more  slowly  with  depth, 
resulting  in  smaller  gradients. 

Another  mechanism  that  would  account  for  these  velocities  involves  the 
reduction  of  porosity.  In  regions  with  a  high  rate  of  deposition  sediments  may 
accumulate  so  rapidly  that  porosity  does  not  have  time  to  be  fully  reduced 
under  overburden  pressure.  The  lowest  values  of  gradients  in  surficial  sediments 
occur  in  such  areas  [22,  pg.  920],  When  porosity  finally  is  reduced,  there  would 
be  a  rapid  increase  in  velocity,  causing  a  larger  gradient.  A  collapse  of  the 
cardhouse  structure  and  the  resulting  rapid  increase  of  intergranular  contacts 
could  also  play  a  role.  After  this  occurs,  gradients  decrease  with  increasing 
depth,  eventually  reaching  a  nominal  value. 

Numerous  reflection  profiles  on  the  Monterey  fan  exist,  but  correlating 
reflection  profiles  with  refraction  profiles  can  be  misleading.  Reflection  data 
indicate  changes  in  acoustic  properties  while  refraction  data  give  velocity  infor¬ 
mation  A  relatively  thin  layer  of  a  different  material  may  give  a  strong 
reflection,  but  would  not  necessarily  affect  the  velocity  gradients  obtained  by 
refraction  data  A  gradual  change  in  velocity  over  a  relatively  large  distance 
will  produce  a  significant  change  in  the  gradient  from  the  refraction  data,  but 
would  show  up  weakly,  if  at  all,  on  a  reflection  profile.  Given  the  complex  evolu¬ 
tion  of  deep  sea  fans,  subsurface  topography  is  bound  to  play  some  role,  but 
without  a  deep  core  from  this  area  it  is  difficult  to  know  which  of  the  above  is 
the  more  likely  explanation. 
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4.  Attenuation  Estimation 


4.1  Introduction 


An  estimate  of  the  attenuation  profile  will  be  made  using  a  method  of 
spectral  ratios  (see  Appendix  A-2).  The  attenuation  coefficient  is  proportional 
to  the  slope  of  the  spectral  ratio  (spectrum  of  the  sediment  refracted  arrival 
divided  by  the  spectrum  of  the  direct  path  arrival).  The  sediment  will  again  be 
modeled  as  a  horizontally  layered  medium,  with  each  layer  having  a  separate 
attenuation  coefficient.  Solving  for  each  of  the  attenuation  coefficients  gives 
attenuation  as  a  function  of  depth. 

An  attenuation  peak  in  the  sediment  has  been  predicted  but  few 
single-experiment  studies  have  been  able  to  confirm  this  peak  19,  29,  30,  45,  63], 
The  results  presented  here  show  an  attenuation  peak  at  around  400-500  meters 
depth  that  is  related  to  a  change  in  sediment  properties  at  that  depth. 

Hamilton  has  compiled  data  from  a  number  of  sources  which  indicate 
an  approximate  linear  dependence  of  the  attenuation  coefficient  k  on  frequency 
19l,  and  that  velocity  dispersion  is  negligible.  Recent  work  has  indicated  that 
this  may  not  be  the  case  [28,  62].  The  issue  of  frequency  dependence  is  not 
addressed  here,  and  it  will  be  assumed  that  k  is  a  linear  function  of 
frequency  (Q-1  is  independent  of  frequency). 

Several  mechanisms  that  can  affect  wave  amplitude  are  discussed  in 
the  next  section. 


4.2  Attenuation  in  marine  sediments 


4-S.l  Factors  affecting  wave  amplitude 
Frictional  losses  [18,  19] 

Frictional  losses  arise  because  of  movement  between  the  sediment  par¬ 
ticles.  Internal  friction  varies  with  grain  size,  the  number  and  kind  of  grain  con¬ 
tacts,  cohesion,  and  friction  between  silt  and  clay  particles.  The  loss  for  this 
type  of  attenuation  is  linearly  related  to  frequency.  It  is  considered  to  be  the 
dominant  loss  mechanism  for  sound  propagation  through  the  sediment.  Porosity 
correlates  well  with  attenuation  in  sediments  since  it  is  a  good  measure  of  the 
number  of  interparticle  contacts. 


Viscoelastic  losses  611 

Viscoelastic  losses  occur  when  the  pore  water  moves  with  respect  to  the 
mineral  frame.  Permeability  of  the  sediment  and  viscosity  of  the  interstitial 
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water  are  important  factors  for  this  type  of  attenuation.  Attenuation  is  propor¬ 
tional  to  f2  for  low  frequencies  (less  than  10  Hz)  and  V/  for  high  frequencies 
(greater  than  lOOKHz). 


Scattering  [61] 

Rayleigh  scattering  occurs  as  wavelengths  approach  grain  size. 
Attenuation  due  to  Rayleigh  scattering  is  proportional  to  /4.  This  should  not  be 
be  a  factor  in  refraction  studies  because  of  the  long  wavelengths  associated  with 
the  low  frequencies  that  are  used  (at  100  Hz,  the  wavelength  is  about  15 
meters). 


Geometrical  spreading  [8,  50] 

Spreading  losses  are  not  losses  that  result  in  an  irreversible  of  conver¬ 
sion  to  heat,  but  cause  wave  amplitude  to  decrease  by  redistributing  the  energy 
over  a  larger  area  In  general,  a  spherical  spreading  loss  model  does  not  apply 
in  sediments  since  the  velocity  changes  as  a  function  of  depth  (although  in  prin¬ 
ciple,  if  the  velocity  profile  is  known  the  spreading  loss  can  be  calculated). 
Losses  due  to  geometrical  spreading  are  not  frequency  dependent. 


Reflection  or  transmission  coefficients  [50] 

Layering  of  the  sediment  will  cause  a  fraction  of  the  energy  to  be 
reflected.  The  more  energy  that  is  reflected,  the  smaller  the  amplitude  of  the 
waveform  becomes.  Thus,  the  magnitude  of  reflection  or  transmission 
coefficients  has  an  effect  on  the  amplitude  of  the  arrival.  This  effect  is  indepen¬ 
dent  of  frequency,  but  may  be  dependent  on  the  angle  of  incidence  as  it  crosses 
each  layer.  It  is  not  always  clear  whether  the  sediment  should  be  interpreted  as 
consisting  of  many  fine  layers,  a  few  thick  layers,  or  something  in  between. 
Energy  reflected  at  an  interface  is  not  lost  through  conversion  to  heat,  and  is 
not  frequency  dependent,  but  in  the  presence  of  many  layers  there  may  appear 
to  be  a  frequency  dependent  loss  due  to  the  effect  of  intrabed  multiples. 


Intrabed  multiples  [50,  54,  55,  56,  58,  59] 

Intrabed  multiples  refer  to  multiple  reflections  off  layers  that  occur  as 
the  wavefront  propagates  through  the  medium.  When  an  initially  downgoing 
pulse  encounters  a  change  in  acoustic  impedance,  part  of  its  energy  is  reflected 
back  off  of  the  interface.  This  reflected  (upgoing)  energy  encounters  a  previous 
layer,  and  a  portion  of  it  is  reflected  back  downward.  Since  it  has  now  been 
reflected  twice,  its  phase  is  the  same  as  that  of  the  original  downgoing  wave.  In 
the  presence  of  many  layers,  each  of  the  multiple  reflections  will  arrive  slightly 
later  in  time  than  the  original  pulse.  This  has  the  effect  of  broadening  the 
received  arrival.  The  more  layering,  the  more  the  energy  is  delayed,  and  the 
more  the  pulse  is  spread  out  in  time.  This  smoothing  and  broadening  of  the 
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pulse  tends  to  tilt  the  spectrum,  shifting  energy  at  higher  frequencies  to  the 
lower  end  of  the  spectrum.  In  this  sense,  the  effect  of  intrabed  multiples  appears 
to  be  frequency  selective.  While  there  is  no  energy  loss  from  conversion  to  heat, 
the  apparent  attenuation  from  intrabed  multiples  is  indistinguishable  from  other 
loss  mechanisms,  such  as  frictional  or  viscous  losses,  which  attenuate  higher  fre¬ 
quencies  more  than  lower  ones.  O’Doherty  and  Anstey  [50],  and  Schoenberger 
and  Levin  [55,  56  show  that  intrabed  multiples  may  account  for  as  much  as  1/3 
to  1/2  of  the  total  attenuation.  Other  references  [54,  58,  59]  also  discuss  this 
important  phenomenon. 


4-2.2  Intrinsic  and  effective  attenuation 

The  previous  section  underscores  the  need  for  defining  just  what  is 
being  measured.  Attenuation  is  often  defined  as  the  decay  in  wave  amplitude 
caused  by  friction  that  results  in  an  irreversible  conversion  of  energy  to  heat. 
This  is  also  referred  to  as  intrinsic  attenuation.  Other  effects  (scattering, 
intrabed  multiples),  which  result  in  a  re-distribution  rather  than  a  loss  of  energy 
also  cause  a  decrease  in  wave  amplitude.  This  type  of  loss  is  called  "apparent 
attenuation.”  The  loss  due  to  both  of  these  effects  is  called  "effective  attenua¬ 
tion.  Although  intrinsic  attenuation  is  a  more  fundamental  parameter,  it  is 
effective  attenuation  that  is  used  in  modeling  sound  propagation. 

The  combination  of  effects  that  make  up  effective  attenuation  compli¬ 
cates  the  spectrum.  Intrinsic  attenuation  is  primarily  due  to  frictional  losses 
and  has  a  fl  dependence.  The  effect  of  intrabed  multiples  is  to  tilt  the  spec¬ 
trum,  raising  the  lower  end  and  decreasing  the  upper  end.  This  tilt  is  not  neces¬ 
sarily  linear  with  frequency,  so  the  spectrum  of  the  combination  may  or  may  not 
show  a  linear  dependence  on  frequency.  Separation  of  intrinsic  from  effective 
attenuation  is  a  research  topic  currently  receiving  attention  [40,  43,  44], 


4-2.3  Previous  studies 

1  reviously  reported  attenuation  values  [19  have  been  obtained  from  a 
variety  of  experiments.  These  include  laboratory  studies,  in  situ  measurements 
made  using  probes  in  the  upper  two  meters  of  sediment,  and  reflection  and 
refraction  experiments.  Results  from  laboratory  experiments  give  attenuation 
for  frequencies  in  the  KHz  to  MHz  range.  At  the  velocities  of  natural  marine 
sediments  (greater  than  1500  m/s)  wavelengths  would  be  on  the  order  of  15 
meters  at  100  Hz  -  very  long  for  a  laboratory  experiment.  Since  small,  homo¬ 
geneous  samples  of  sediment  are  used,  these  studies  give  a  measure  of  intrinsic 
attenuation.  Attenuation  at  high  frequencies  is  found  from  reflection  studies, 
but  depths  to  as  much  as  a  kilometer  are  sampled,  so  effective  attenuation  is 
found.  Refraction  studies  provide  information  about  attenuation  at  low  fre¬ 
quency.  Large  distances  and  vertically  inhomogeneous  sections  of  sediment  are 
covered,  so  these  studies  measure  effective  attenuation.  A  summary  of  the  types 
of  studies  is  given  in  Table  4.1. 
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Attenuation  Studies 

Type  of  Study 

Frequencies 

Used1 

Type  of  Atten¬ 
uation  Measured 

2 

Region  of  Sediment^ 

in  situ  (probes) 

high  frequency 

intrinsic 

surface 

reflection 

high  frequency 

effective 

shallow  to  moderately  deep 

refraction 

low  frequency 

effective 

moderately  deep  to  deep 

1  low  frequency  <300  Hz,  high  frequency  >1  kHz 

2  shallow  <2  meters,  deep  >500  meters 


Table  4*1  Attenuation  studies. 


Studies  of  attenuation  generally  yield  attenuation  as  a  function  of  fre¬ 
quency,  or  attenuation  as  a  function  of  depth 


4.2.4  Attenuation  as  a  function  of  frequency 

Frictional  losses  are  thought  to  be  the  dominant  attenuative  mechan- 
lsm-  Hamilton  |19  has  compiled  data  from  a  number  of  sources  and  concludes 
that  k  is  approximately  linearly  dependent  on  frequency.  An  important  factor 
which  would  result  in  a  nonlinear  dependence  on  frequency  is  the  movement  of 
pore  water  with  respect  to  the  mineral  frame.  Movement  of  the  pore  water 
results  in  viscous  losses  and  velocity  dispersion.  Velocity  dispersion  would  be 
largest  in  sands  since  permeability  is  high  and  there  would  be  a  higher  amount 
of  movement  of  pore  water.  In  reviewing  studies  of  dispersion  in  sands,  Hamil¬ 
ton  concluded  that  dispersion  is  negligible  from  a  few  KHz  to  1  MHz  in  sands 
and  less  than  1  KHz  to  2  MHz  in  silt-clays,  although  he  points  out  that  the  tests 
were  made  over  an  order  of  magnitude  of  frequency  or  less,  which  may  not  be 
enough  of  a  range  to  show  dispersion.  Recent  reports  suggest  that  a  strict  linear 
dependence  on  frequency  is  not  correct  [28,  62], 


4-2.5  Attenuation  as  a  function  of  depth 

There  is  little  experimental  data  on  attenuation  as  a  function  of  depth 
in  the  sediments.  Hamilton  [19]  predicted  an  initial  increase  in  attenuation  to 
around  200  meters  in  depth  followed  by  a  decrease  in  attenuation  with  increas¬ 
ing  depth  His  prediction  is  based  on  two  opposing  effects:  (l)  Attenuation 
increases  as  porosity  decreases.  Porosity  is  defined  as  the  volume  of  water-filled 
pore  space  per  volume  of  sediment.  As  porosity  decreases,  sediment  grains  come 
into  closer  contact,  which  increases  the  number  of  intergranular  contacts  and 
the  friction  between  them.  This  increase  in  friction  causes  an  increase  in 
attenuation.  (2)  Attenuation  is  decreased  as  overburden  pressure  increases.  As 
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depth  (pressure)  increases,  the  sediment  grains  are  forced  together  and  begin  to 
form  a  more  rigid  structure.  There  is  less  movement  of  the  grains  and  thus  less 
frictional  loss.  Attenuation  therefore  decreases  with  increasing  overburden  pres¬ 
sure. 

Silt-clays  are  the  dominant  part  of  surficial  sediments  on  deep  sea  fans. 
The  porosity  of  these  sediments  is  between  70  and  90*^5.  There  is  a  large  initial 
decrease  in  porosity  with  increasing  depth,  and  that  effect  dominates,  causing  an 
increase  in  attenuation.  By  100  to  200  meters,  there  is  little  further  reduction  in 
porosity  and  the  effect  of  increasing  overburden  pressure  begins  to  dominate. 
Past  this  point,  there  is  a  decrease  in  attenuation  with  increasing  depth.  Hamil¬ 
ton  reports  a  value  for  k  at  the  surface  between  .05  to  .1  dB/(km-Hz),  followed 
by  an  increase  to  .2  dB  /(km-  Hz)  at  100-200  meters,  and  a  decrease  in  k  past  200 
meters.  (Attenuation  coefficient  k  is  defined  by  equation  A2,5a  or  A2.5c  in 
Appendix  A2.  It  has  units  of  dB /(meter- kHz),  or  equivalently,  dB /(km-Hz  ) 

Jacobson,  Shor,  and  Dorman  [30]  offer  a  different  explanation  for  the 
shape  of  the  attenuation  curve.  They  found  a  peak  of  attenuation  at  600 
meters,  and  feel  that  the  observed  peak  is  a  consequence  of  sediment 
hthification.  600  meters  is  approximately  the  depth  at  which  the  unlithified  sed¬ 
iment  becomes  a  mudstone.  The  collapse  of  the  cardhouse  structure  occurring 
at  that  depth  increases  intergranular  contacts  and  aids  in  cementation  of  the 
mineral  grains.  This  reduces  sliding  friction,  which  is  the  dominant  mechanism 
for  attenuation  in  the  sediment. 

They  also  point  out  that  the  decrease  of  k  with  increasing  depth  that 
Hamilton  reports  is  more  likely  due  to  the  units  of  k  than  to  overburden  pres¬ 
sure.  The  attenuation  coefficient  k  gives  a  measure  of  the  loss  per  distance.  As 
the  velocity  increases  with  depth,  there  is  a  concomitant  increase  in  the 
wavelength  The  pathlength  over  one  cycle  therefore  becomes  longer  and  longer. 
A  unit  of  length  becomes  a  smaller  and  smaller  portion  of  a  wavelength,  and  k 
decreases.  The  quality  factor  Q"1  measures  attenuation  as  loss  per  cycle'  and  is 
thus  independent  of  frequency  and  velocity.  For  this  reason,  they  feel  Q-1  is 
preferable  to  k  as  a  measure  of  attenuation. 


4-2.6  Interpretation  of  results 

The  data  that  have  been  compiled  appear  to  have  a  somewhat  remark¬ 
able  linear  relationship  over  a  wide  range  of  frequencies  [19,  Fig.  1],  These 
results  were  obtained  from  different  types  of  experiments  as  well  as  different  geo¬ 
graphical  areas.  A  compilation  of  data  from  different  regions  and  sediment 
types,  and  obtained  using  different  techniques,  must  contain  a  large  amount  of 
scatter.  Scatter  results  from  the  estimation  of  different  quantities  (intrinsic  or 
effective),  and  because  attenuation'  is  depth  dependent.  Given  the  wide  variety 
of  sources  for  data,  and  the  amount  of  scatter,  inferences  drawn  from  them  need 
to  be  interpreted  carefully. 
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4.3  Methods  of  estimating  attenuation 


Several  methods  have  been  used  to  estimate  attenuation  in  a  medium. 
These  can  be  grouped  into  time  domain  and  frequency  domain  methods.  Fre¬ 
quency  domain  methods  consider  some  aspect  of  the  source  and  attenuated  spec¬ 
tra.  Time  domain  methods  use  a  measure,  such  as  the  change  in  the  peak 
amplitude  of  the  waveform,  to  relate  changes  in  the  time  series  with  attenua¬ 
tion.  Without  considering  the  frequency  domain  characteristics  of  the  signals 
(source  and  received  waveforms),  time  domain  methods  generally  produce 
misleading  or  incorrect  results  An  example  will  show  why  this  is  so 

Consider  a  refraction  experiment  that  uses  two  different  sources.  If  the 
first  source  is  broadband  (an  input  spike,  which  has  a  flat  spectrum,  or  an  explo¬ 
sion,  which  has  a  relatively  broadband  spectrum),  much  of  the  energy  in  the  sig¬ 
nal  is  carried  in  the  higher  frequencies.  Since  high  frequencies  are  highly 
attenuated,  much  of  the  high  frequency  energy  will  be  lost  and  the  waveform 
will  appear  to  decay  fairly  rapidly.  If  a  second  source  which  produces  a  low  fre¬ 
quency  narrow-band  signal  is  then  used,  the  energy  that  was  previously 
transmitted  at  higher  frequencies  is  now  carried  at  lower  frequencies.  In  this 
case,  since  there  is  relatively  low  attenuation  at  the  lower  frequencies,  the 
waveform  will  show  less  attenuation  than  when  the  broadband  source  was  used. 
Although  the  acoustic  energy  travels  the  same  distance  through  the  same 
medium,  the  amplitude  of  the  received  waveform  is  different  depending  upon  the 
source.  This  example  shows  that  the  frequency  characteristics  of  the  source 
must  be  considered  in  these  types  of  experiments. 

A  method  of  spectral  ratios  is  used  in  the  present  study.  Spectral 
ratios  are  formed  by  dividing  the  spectrum  of  the  refracted  path  by  that  of  the 
direct  path,  thus  taking  the  source  (direct  path)  spectrum  into  account. 
Attenuation  is  proportional  to  the  slope  of  the  spectral  ratio.  The  method  is 
derived  in  Appendix  A2. 


4.4  Data 

Figures  4.1a  through  4.1d  give  examples  of  the  spectral  ratio  data.  The 
top  plots  in  these  figures  are  the  hydrophone  time  series,  the  middle  plots  show 
the  spectra  of  the  two  time  series,  and  the  bottom  plots  give  the  spectral  ratios. 
A  linear  least  squares  fit  to  each  spectral  ratio  is  also  shown.  Spectral  ratio 
slopes,  needed  in  equation  A2.13  (Appendix  A2),  are  obtained  from  the  least 
squares  fits. 

The  effect  of  increasing  path  length  can  be  seen  in  these  figures.  Each 
figure  in  the  sequence  presents  data  from  a  deeper  ray,  or  equivalently,  from  a 
ray  which  has  a  longer  path  length.  Arrivals  which  have  longer  path  lengths 
suffer  more  attenuation,  and  the  refracted  spectra  have  substantially  lower  lev¬ 
els  than  their  respective  source  spectra.  The  result  is  that  spectral  ratio  slopes 
are  steeper  for  arrivals  that  penetrate  deeper  in  the  sediment. 


58 


SIO  Reference  87-28 


D  Shot  41  BR 
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Figure  4.1a  Spectral  ratio  data.  Shot  41  BR. 

Time  series  at  top  are  scaled  to  have  same  maximum  amplitude. 
Data  is  from  a  shallow  refraction  (100  meters  depth  in  sediment). 
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D  Shot  17 


BR 


frequency  (Hz) 


slope  =  -.042 


dB 

Hz 


Figure  4.1b  Spectral  ratio  data.  Shot  17  BR. 

Time  series  at  top  are  scaled  to  have  same  maximum  amplitude. 

Data  is  from  a  moderately  shallow  refraction  (133  meters  depth  in  sediment). 
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Figure  4.1c  Spectral  ratio  data.  Shot  29  S-BR-S. 

Time  series  at  top  are  scaled  to  have  same  maximum  amplitude. 

Data  is  from  a  medium-depth  refraction  (213  meters  depth  in  sediment). 
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D  Shot  44  S-BR-S 


frequency  (Hz) 


slope  =  -.190  — 
Hz 


Figure  4.1d  Spectral  ratio  data.  Shot  44  S-BR-S. 

Time  series  at  top  are  scaled  to  have  same  maximum  amplitude. 
Data  is  from  a  deep  refraction  (470  meters  depth  in  sediment). 
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Spectral  ratio  data 

Floor 

Depth  in 

Path 

Total 

-slope 

-slope 

Arrival 

Depth 

Sediment 

Length 

Time 

of  SRt 

of  SR0 

17a  br 

2822 

133 

1510 

.9760 

.042194 

.004858 

19a  br 

2830 

88 

1216 

.7971 

.029922 

.003445 

19b  brs 

2835 

455 

2876 

1.6553 

.156992 

.018074 

21b  brs 

2850 

246 

1912 

1.1887 

.076800 

.008842 

23b  brs 

2860 

191 

1682 

1.0667 

.028160 

.003242 

23b  sbrs 

2833 

516 

3186 

1.7896 

.203408 

.023418 

25a  sbr 

2820 

165 

1578 

1.0106 

.099072 

.011406 

25a  brs 

2875 

145 

1547 

.9970 

.065424 

.007532 

25b  sbrs 

2843 

323 

2196 

1.3279 

.130222 

.014992 

27a  sbr 

2823 

126 

1466 

.9498 

.056320 

.006484 

27a  brs 

2887 

115 

1394 

.9066 

.044800 

.005158 

29a  sbrs 

2856 

213 

1772 

1.1152 

.091622 

.010548 

31a  sbrs 

2862 

186 

1661 

1.0554 

.054918 

.006323 

33a  sbrs 

2865 

169 

1595 

1.0195 

.055712 

.006414 

36a  br 

2848 

47 

877 

.5829 

.026532 

.003055 

36a  sbr 

2817 

202 

1725 

1.0903 

.026338 

.003032 

36a  brs 

2863 

181 

1641 

1.0445 

.033556 

.003863 

36a  sbrs 

2835 

470 

2951 

1.6904 

.163328 

.018804 

40a  br 

2821 

154 

1546 

.9939 

.024320 

.002800 

41a  br 

2827 

100 

1297 

.8471 

.020000 

.002303 

41a  brs 

2830 

622 

3611 

1.9807 

.196820 

.022660 

43a  sbr 

2814 

272 

2014 

1.2406 

.040872 

.004706 

43a  brs 

2851 

233 

1859 

1.1612 

.092058 

.010599 

44a  sbr 

2816 

207 

1748 

1.1026 

.045708 

.005262 

44a  brs 

2861 

184 

1651 

1.0499 

.055792 

.006423 

44a  sbrs 

2835 

470 

2951 

1.6904 

.190280 

.021907 

45a  br 

2850 

43 

842 

.5605 

.018530 

.002133 

45a  sbr 

2819 

174 

1612 

1.0290 

.035788 

.004120 

45a  brs 

2868 

158 

1557 

.9994 

.039040 

.00495 

46a  brs 

2875 

141 

1554 

1.0024 

.025600 

.002947 

47a  sbr 

2823 

122 

1437 

.9324 

.027686 

.003187 

47a  brs 

2885 

115 

1394 

.9066 

.018022 

.002075 

47a  sbrs 

2850 

265 

1995 

1.2315 

.076800 

.008842 

48a  sbrs 

2855 

218 

1796 

1.1280 

.051200 

.005895 

Table  4.2  Spectral  ratio  data.  Floor  depth  is  the  depth  at  which  the  ray  in¬ 
tersects  the  sea  floor.  Path  length  is  total  path  length  in  the  sediment.  SRk 
refers  to  the  definition  of  the  spectral  ratio  given  in  equation  A2.6.  SRQ  refers 
to  the  definition  of  the  spectral  ratio  given  in  equation  A2.9. 
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Figure  4  la  is  data  from  a  ray  which  penetrated  to  a  depth  of  100 
meters  and  had  a  path  length  of  1297  meters.  In  Figure  4-lb  the  ray  turned  at 
at  depth  of  133  meters  and  had  a  path  length  of  1510  meters.  The  increase  in 
path  length  of  213  meters  has  doubled  the  spectral  ratio  slope.  Figure  4. Id  is 
data  from  a  ray  which  had  a  depth  of  470  meters  and  a  path  length  of  2951 
meters.  The  spectral  ratio  slope  is  nearly  ten  times  that  of  the  slope  in  Figure 
4  la.  In  general,  spectral  ratios  of  shallow  refractions  have  a  linear  region  at 
higher  frequencies,  and  are  linear  over  a  larger  frequency  interval,  than  highly 
attenuated  deeper  refractions.  Deep  refractions  lose  most  of  their  high  frequency 
energy  and  show  linearity  over  smaller  intervals  since  the  signal  at  high  frequen¬ 
cies  is  down  in  the  noise  floor,  and  the  spectrum  flattens  out.  Spectral  ratio 
slopes  (557?)  along  with  the  total  time  and  path  length  for  each  ray  are  given  in 
Table  4.2. 

The  spectra  used  to  calculate  the  spectral  ratios  are  found  by  separat¬ 
ing  individual  arrivals  from  the  time  series  using  a  rectangular  window.  A  fast 
Fourier  transform  (FFT)  of  each  channel  for  a  particular  ray  type  (i.e.  D  BR, 
etc.)  is  calculated,  and  the  spectra  are  incoherently  averaged  (magnitude 
squared).  The  square  root  of  the  average  is  then  taken  to  obtain  the  magnitude 
spectrum. 

The  seven  good  channels  from  the  low  gam  half  of  the  array  are  used 
when  possible.  In  many  instances  the  BR-S  and  S-BR  arrivals  cross  in  time  and 
only  two  or  three  of  the  channels  can  be  used.  For  example,  in  Figure  2.5c, 
these  two  arrivals  converge  near  hydrophone  II  There  is  complete  overlap  by 
hydrophone  8  (not  shown  in  this  figure)  and  the  arrivals  separate  again  lower  in 
the  array.  In  these  cases,  only  channels  in  which  arrivals  are  separated  in  time 
are  used  As  shot  range  increases,  the  direct  (D)  and  BR  paths  run  together. 
When  these  two  arrivals  begin  to  overlap,  the  surface  reflected  path  (S)  is  used 
to  find  the  source  spectrum. 


Solutions  presented  in  Fieure  4.2 

Figure 

number 
of  layers 

eigenvectors 

retained 

SSR 

X  10-4 

4.2a 

3 

3 

1  85 

4.2b 

4 

4 

1.61 

4.2c 

7 

5 

1.53 

4. 2d 

13 

5 

2,09 

4.2e 

13 

7 

1.60 

Table  4,3  Solutions  presented  in  Figure  4.2. 
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Figure  4.2a  Attenuation  profile:  3  layer  model,  3  eigenvectors  retained. 

A  theoretical  curve  generated  from  the  attenuation  profile  is 
plotted  along  with  the  experimental  data. 
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Figure  4.2b  Attenuation  profile:  4  layer  model,  4  eigenvectors  retained. 

A  theoretical  curve  generated  from  the  attenuation  profile  is 
plotted  along  with  the  experimental  data. 
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Figure  4.2c  Attenuation  profile:  7  layer  model,  5  eigenvectors  retained. 

A  theoretical  curve  generated  from  the  attenuation  profile  is 
plotted  along  with  the  experimental  data. 
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Figure  4. 2d  Attenuation  profile:  13  layer  model,  5  eigenvectors  retained. 

A  theoretical  curve  generated  from  the  attenuation  profile  is 
plotted  along  with  the  experimental  data. 
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Figure  4.2 e  Attenuation  profile:  13  layer  model,  7  eigenvectors  retained 
A  theoretical  curve  generated  from  the  attenuation  profile  is 
plotted  along  with  the  experimental  data. 
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4.5  Results 

The  spectral  ratio  slopes  from  Table  4.2  were  used  in  equation  A2.13 
and  inverted  using  a  singular  value  decomposition  (SVD)  [15,  39,  41,  66],  Four 
models,  with  3,  4,  7,  and  13  layers  are  presented,  and  are  summarized  in  Table 
4.3.  Results  are  illustrated  in  Figures  4,2a  through  4.2e.  The  top  plot  in  these 
figures  gives  Q  as  a  function  of  depth  in  the  sediment.  In  the  bottom  plot,  a 
theoretical  curve  (solid  line)  generated  from  equation  A2.15  is  compared  with 
the  experimental  data. 

The  3  and  4  layer  models  had  well  conditioned  matrices  and  all  eigen¬ 
vectors  were  incorporated  into  the  solution.  Solutions  became  unstable,  produc¬ 
ing  negative  attenuation  coefficients,  when  all  eigenvectors  for  the  7  layer  model 
and  more  than  8  eigenvectors  for  the  13  layer  model  were  included  in  the  solu¬ 
tion  The  sum  of  squared  residuals  (  SSR  =  £(<557?, --^SR,.)2 )  gives  an  indication 

of  the  quality  of  fit.  The  7  layer  model  (Figure  4  2c)  has  the  lowest  SSR  (last 
column  in  Table  4.3),  suggesting  that  we  are  able  to  obtain  better  resolution 
than  is  given  by  a  3  or  4  layer  model.  Beyond  this,  however,  the  values  increase 
because  of  data  resolution  problems.  This  effect  is  seen  in  the  13  layer  model. 
Retaining  5  eigenvectors  produces  a  smooth  attenuation  profile,  but  does  not 
adequately  represent  the  data  (Figure  4.2d).  Retaining  8  eigenvectors  results  in 
a  better  fit  to  the  experimental  data,  but  produces  an  attenuation  profile  with  a 
large  amount  of  variance  (Figure  4.2e). 

In  selecting  a  solution,  we  can  reject  the  3  layer  model  in  favor  of  the  4 
layer  model  based  on  an  F-test.  Although  the  sum  of  squared  residuals  (SSR)  is 
less  for  the  7  layer  model  than  for  the  4  layer  model,  there  are  fewer  degrees  of 
freedom,  the  reduction  in  b!SR  is  not  large  enough  to  warrant  rejection  of  the  4 
layer  model  using  an  F-test.  The  4  layer  model  along  with  90%  confidence 
bounds  (bold  lines)  is  shown  in  Figure  4.3 


Q'1 

Figure  4.3  Four  layer  attenuation  model  with  90%-  confidence  bounds. 

Model  is  given  by  the  light  line,  confidence  bounds  are  given  by  bold  lines. 
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Although  the  4  layer  model  cannot  be  rejected  in  favor  of  the  7  layer 
model,  this  does  not  mean  that  the  7  layer  model  must  be  rejected.  It  can  be 
shown  (also  by  use  of  an  F-test)  that  the  variances  of  these  two  models  are  not 
significantly  different.  The  7  layer  model  is  somewhat  more  appealing,  as  it  pro¬ 
vides  better  resolution  while  still  yielding  a  smooth  solution.  For  these  reasons, 
both  the  4  and  7  layer  solutions  are  considered  acceptable.  Details  of  these  two 
solutions  are  given  in  Tables  4  4  and  4.5. 


Attenuation  profile  -  4  layer  model 

Depth  (meters) 

k  (dB/ m/kHz) 

Q (nepers/radian) 

0-200 

.027 

.0015 

200-400 

.057 

.0036 

400-600 

.081 

.0057 

600-800 

.018 

.0015 

Table  4  4  Attenuation  profile  for  a  4  layer  model 
Attenuation  profile  shown  in  Figure  4.2b  (4  eigenvectors  retained). 


Attenuation  profile  -  7  layer  model 

Depth  (meters) 

k  (dB/ m/kHz) 

Q-1  (nepers/radian) 

0-100 

.022 

.0012 

100-200 

.033 

.0019 

200-300 

.052 

.0032 

300-400 

.057 

.0040 

400-500 

.085 

.0057 

500-600 

.071 

.0052 

600-700 

.040 

.0028 

Table  4.5  Attenuation  profile  for  a  7  layer  model 
Attenuation  profile  shown  in  Figure  4.2c  (5  eigen- 
ectors  retained). 


4.6  Discussion 

The  most  notable  feature  in  all  solutions  is  an  attenuation  peak  at 
about  400-500  meters  depth.  Hamilton  predicted  an  attenuation  peak  due  to 
the  opposing  effects  of  porosity  reduction  and  the  increase  in  overburden 
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pressure  as  depth  increases,  while  Jacobson,  Shor,  and  Dorman  [30],  and  Jacob¬ 
son,  Shor,  and  Bee  [29]  felt  that  the  peak  they  found  was  a  result  of  sediment 
lithification  (see  the  discussion  in  Section  4.2.5).  Without  a  deep  core,  it  is 
difficult  to  know  which  of  the  two  effects  is  responsible. 

The  relationship  between  the  attenuation  profile  and  the  velocity  gra¬ 
dient  is  shown  in  Figure  4.4  A  four  layer  velocity  profile  for  the  sediment  was 
presented  in  Table  3  8.  The  gradients  from  this  result  are  plotted  along  with 
the  attenuation  profile  for  a  four  layer  model.  In  this  figure,  it  is  seen  that  there 
are  peaks  in  both  the  velocity  gradient  and  attenuation  profiles.  Possible  causes 
for  a  peak  in  the  velocity  gradient  have  been  discussed  in  Section  3  8.  The  velo¬ 
city  gradient  and  attenuation  are  affected  similarly  by  physical  changes  occur¬ 
ring  in  the  sediment;  the  discussion  in  that  section  also  applies  to  the  attenua¬ 
tion  profile.  A  difference  between  these  two  profiles  is  that  the  attenuation 
peaks  at  a  greater  depth  than  the  velocity  gradient.  The  attenuation  peak  also 
occurs  below  the  velocity  gradient  peak  in  Jacobson,  Shor,  and  Dorman  [30, 
Figs.  4  and  9],  and  in  Jacobson,  Shor,  and  Bee  [  29,  Figs.  9  and  16  . 

Changes  in  velocity  and  attenuation  in  the  upper  sediments  (before 
lithification  begins)  are  primarily  due  to  changes  in  porosity  [19,  22],  The  high 
porosity  (70-90%)  silts  and  clays  present  at  the  surface  of  deep  sea  fans  are 
rapidly  reduced  with  increasing  depth  in  the  sediment.  The  largest  amount  of 
reduction  in  porosity  is  thought  to  occur  by  200-300  meters  [19,  20 

Velocity  increases  monotonically  as  porosity  decreases  [25,  Fig.  6;  As 
the  rate  of  change  in  porosity  begins  to  slow,  the  velocity  gradients  would  begin 
to  approach  a  constant.  In  the  data  presented  in  Figure  4.4  this  occurs  around 
400  meters. 

Attenuation  is  not  a  monotonic  function  of  porosity  [19,  Fig.  2  ]. 
Attenuation  increases  as  porosity  decreases,  reaching  a  peak  between  50-55%. 
As  porosity  decreases  below  50%,  attenuation  decreases.  If  the  porosity  were 
still  greater  than  55  or  60  %  at  400  meters,  then  a  further  reduction  in  porosity 
(as  depth  increases)  would  result  in  an  increase  in  attenuation.  If  changes  in 
porosity  were  the  only  effect,  then  a  peak  in  the  velocity  gradient  gives  an  indi¬ 
cation  of  where  the  most  rapid  reduction  in  porosity  occurs;  the  attenuation 
peak  would  indicate  the  depth  where  porosity  values  are  around  50%  (or  the 
place  where  effects  of  overburden  pressure  counterbalance  the  reduction  in 
porosity).  Other  factors  are  most  likely  present,  with  lithification,  no  doubt, 
playing  an  important  role.  As  sediments  lithify,  the  mineral  frame  becomes 
rigid,  and  the  attenuation  and  velocity  gradient  will  approach  a  constant  value. 

Finally,  in  Figures  4.5  and  4.6,  the  results  of  this  study  are  added  to 
the  attenuation  data  compiled  by  Hamilton  [19].  Attenuation  data  as  a  func¬ 
tion  of  frequency  in  Figure  4.5  is  plotted  as  a  band  of  values.  This  band 
corresponds  to  the  smallest  and  largest  values  found  (Table  4.5).  Attenuation  as 
a  function  of  depth  is  shown  in  Figure  4.6. 


72 


700  600  500  400  300  200  100 


SIO  Reference  87-28 


velocity  gradient  (sec'1) 

■40  .55  JO  .85  1.00  1.15 


Q1 


Figure  4.4  Relationship  between  the  velocity  gradient  and  1/Q  as  a  function  of  depth. 
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Figure  4.5  Attenuation  as  a  function  of  frequency. 

Shaded  rectangle  shows  results  of  this  study. 
(Reprinted  from  reference  [19]) 
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Figure  4.6  Attenuation  as  a  function  of  depth. 

(Reprinted  from  reference  [19]) 
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Appendix  Al 


Inverting  X(p),  T(p),  and  r(p)  data  to  obtain  velocity 
as  a  function  of  depth 


Introduction 

In  this  appendix  equations  for  inferring  the  velocity  structure  of  a 
medium  from  T{p\  X(p),  or  r(p)  data  are  developed.  The  end  result  will  be  a 
system  of  linear  equations  which,  when  inverted,  will  yield  information  about  the 
velocity  as  a  function  of  depth. 


A  general  linear  inverse  problem 

A  general  way  of  posing  a  linear  inverse  problem  [64]  is 
V 

d{  =  J  K{(u)  m(u)  du  2=4,2,. ..,7V 

*o 


(All) 


d{  is  the  data  point,  K{  is  the  data  kernel  for  the  ith  datum,  and  m(u)  is  the 
model.  The  kernel  represents  what  is  known  from  the  physics  of  the  problem. 
Given  a  set  of  data  rf,-,  the  goal  is  to  find  the  model  m[u)  which  generates  that 
data. 

One  way  of  proceeding  is  to  define  the  model  as  a  linear  combination 

M 

m(«)  =  E  ai  FAU)  (A1.2) 

where  the  ay  are  unknown  parameters,  and  the  F}(u)  are  specified  in  such  a  way 
as  to  construct  the  desired  model.  Substituting  this  expression  for  the  model 
into  the  general  form  (Al  l)  gives 


M 


Ki{u)  £  a,.  Fy(u) 

y-4 


du 


M  '1 

=  E  ai  I  Ki(u )  FAV )  du  *'=1,2,... ,7V 

>= 1  *o 


(A1.3) 


The  kernel  K{  is  known  and  the  F}-  have  been  specified,  so  the  integral  may  be 
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evaluated  (either  analytically  or  numerically)  for  each  i  and  j.  Let 

•f 


9ij  =  f  K{{ u)  Fj{u)  du 


then  (Al  3)  may  be  written  as 

M 


d\  S  aj  Si]  i—l,2 


(A1.4) 


which  represents  a  set  of  N  linear  equations  in  M  unknowns.  In  matrix  notation 
this  is  d  =Ga.  There  are  N  data  and  M  unknown  model  parameters,  d  is  the 
data  vector,  G  is  a  known  matrix,  and  a  is  the  unknown  vector  of  model 
parameters.  There  are  usually  many  more  data  than  model  parameters,  produc¬ 
ing  a  non-square,  overdetermined  system  of  linear  equations.  An  extensive 
literature  exists  for  solving  this  type  of  problem  [1  (Sect  12.3),  15,  39,  41,  57,  64, 


Ray  propagation  as  an  inverse  problem 

Consider  a  ray  propagating  through  a  medium  whose  velocity  structure 
is  laterally  homogeneous,  but  vertically  inhomogeneous  (figure  Al.l).  11(2),  the 
reciprocal  of  velocity,  is  called  slowness:  u{z)  =  l/v{z).  If  the  velocity  increases 
(slowness  decreases)  with  depth,  then  a  ray  leaving  the  surface  2  =  0  at  angle  d0 
will  refract  as  it  travels  through  the  layers,  returning  to  the  surface  at  a  hor¬ 
izontal  range  of  X  from  where  it  started  The  time  it  takes  to  travel  that  dis¬ 
tance  is  T  seconds.  The  range  X  and  time  T  for  such  a  raypath  is  given  by  [8] 


(A1.5) 


(A1.6) 


Subscripts  on  X  and  T  indicate  the  ith  ray  (f=l,2,...,A).  Subscripts  on  X  and  T 
indicate  the  itk  ray  (*=1,2,...,  IV).  A  commonly  used  linear  combination  of  X(p) 
and  T(p),  called  the  delay  time,  is  defined  by 


t>(p  )  =  T{(p)-pX{(p) 
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f  [  “2(*)  -P2]1/2  dz 


(A1.7) 


In  the  above  equations  p  is  the  ray  parameter,  defined  as 
p  =  u (0)  sin(0o),  where  u(0)  is  the  slowness  at  the  surface  *  =0,  and  60  is  the 
takeoff  angle  (the  angle  the  ray  makes  with  the  vertical  at  the  surface).  If  the 
medium  is  laterally  homogeneous  and  horizontally  layered,  then  p  is  constant 
along  the  entire  raypath:  p  =«  sin(9).  In  this  case  the  ray  parameter  uniquely 
defines  a  ray.  zTP  is  the  depth  at  which  the  ray  turns  (see  figure  Al.l).  Due  to 
symmetry,  the  total  range  and  time  is  twice  the  range  and  time  at  the  turning 
point  (midpoint  of  the  ray  path).  The  limits  of  the  integrals  go  from  z  —  0  to 

*  =  z  TP  >  so  multiplication  by  a  factor  of  two  is  required  to  obtain  values  for  the 
entire  raypath. 

Equations  Al.5,  Al  6,  and  Al.7  define  the  forward  problem.  Given  a 
slowness  profile  u{z)  and  a  ray  parameter  p,  the  time  T ,  horizontal  range  X,  and 
delay  time  r  for  the  raypath  may  be  found  by  evaluating  the  integrals.  '  The 
inverse  problem  is:  given  a  set  of  range  or  time  data,  what  can  be  inferred  about 
the  slowness  (velocity)  profile. 

A  model  which  has  received  considerable  use  is  one  having  constant 
velocity  layers.  A  model  that  has  layers  of  constant  velocity  is  not  very  realistic 
and  often  produces  physically  unreasonable  results.  A  more  realistic  model  may 
be  constructed  if  each  layer  is  assumed  to  have  a  constant  velocity  (slowness) 
gradient  rather  than  a  constant  velocity  (see  figure  Al.2).  Thus,  the  0 /**  layer 

would  have  a  constant  slowness  gradient  of  [du/dz^j.  Use  of  this  model  will  first 
be  considered  for  X(p). 

Before  proceeding,  a  change  of  notation  will  be  made  to  make  the 
discrete  nature  of  the  data  more  explicit.  Given  the  assumption  of  a  laterally 
homogeneous,  horizontally  layered  medium,  each  arrival  is  uniquely  specified  by 
a  ray  parameter  p.  In  view  of  this,  let  the  itk  datum  be  given  by  X(p,).  Then, 
changing  the  variable  of  integration  from  depth  *  to  slowness  «  in  (Al.5)  results 


in 


TP 


(A1.8) 


In  the  notation  of  equation  Al  l  the  kernel  is 


— — '-r—  and  the  model  is  (— ) 
(«  -  Pi)1'2  y  du  ’ 


the  inverse  of  the  slowness  gradient.  The  limits  now  go  from  slowness  at  the 
surface  («0)  to  slowness  at  the  turning  point  (urp)  Next,  as  in  equation  Al.2, 
define  the  model  to  be 


(A1.9) 
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and  uj  =  the  inverse  slowness  gradient  in  the  jtk  layer.  This  essentially 

breaks  the  medium  up  into  M  layers,  modeling  each  layer  with  a  constant,  wy. 
There  are  A/+1  interfaces  (specified  by  u0,  u^,  ...,  u ^  to  define  the  M  layers. 
Slowness  is  assumed  to  decrease  with  depth  (or,  equivalently,  velocity  increases 
with  depth).  Substitution  of  the  model  (Al.9)  into  (Al.8)  yields  a  form  similar 
to  (Al.3) 


X(p, 


■>  -  / 


%Pi 


M 


u2-p?y'2  § F,[u) 


du 


(Al.10) 


The  limits  on  the  integral  extend  only  over  the  jth  layer  since  F}-(u)  =0 
unless  u  e  [ujt  The  integral  may  now  be  evaluated 


du _ 

-a2)1/2 


In  [  u  +  (u2-?,-2)1/2  ] 


-  In  [  Uj  +  («/ -p  2)1/2  ]  -  In  [  +  («/_[  -p2)1/2  ] 

=  in  uy  +  (m/~Pi2)1/2 

11  l  +  KW)1/2  . 


The  limits  are  correct  if  the  ray  crosses  the  jtk  interface.  If  the  ray 
turns  before  the  jtk  interface  (i.e.  if  Pi  >«y)  then  the  upper  limit  should  be  Pi 
instead  of  «y,  otherwise  the  quantity  («?  -  p,?)1/2  becomes  negative.  At  the  turn¬ 
ing  point,  the  ray  is  horizontal  and  9jp  =7 r  .  Since  p  =u  stn(9)  everywhere  along 
the  raypath,  at  the  turning  point  u  =  p  .  Define  «max  =  max(p,-,«y)  so  that  the 
upper  limit  is  either  «y,  the  slowness  at  the  jth  interface  if  the  ray  crosses  it,  or 
p,  the  slowness  at  the  turning  point  if  the  ray  turns  before  the  jth  interface. 
Define  g{]-  to  be  the  quantity  in  brackets  from  equation  Al.10  and  we  have 

*7 

O  f  du 

9<i  ~ 


“max  +  (“L  -P2)1/2 

u]-i  +  (“/- 1  -  P2)1/2 


“max  =  max(p,-,«y)  (Al.ll) 
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Substitution  of  this  result  into  (Al.10)  gives  the  system  of  linear  equations  Al.4. 

x(p<)  =  E  9a 

j=4 

In  matrix  form  this  is 


‘X(p  i)‘ 

Sll  8 12  8 13  •  •  •  8  lAf" 

X(pd 

8 21  822  823  ■  ■  ■  82M 

X(Pz) 

8  31  8  32  8  33  •  •  •  £  2M 

CO] 

0)2 

= 

0)3 

Xips) 

8n  1  #V2  #V3  •  ■  •  8nm 

which  may  be  written  more  compactly  as 

x  =  G  w 


(Al.12) 


Vector  w  contains  the  reciprocals  of  the  unknown  slowness  gradients  (— ).  x 

contains  the  data,  and  matrix  G  is  a  known  matrix  whose  elements  (for  the  itk 
datum  and  the  jtk  layer)  are  given  by  (Al.ll).  The  model  generally  has  a  few 
layers,  while  there  are  many  more  X,  T,  or  t  data.  This  results  in  an  overdeter¬ 
mined  system  of  equations. 

An  expression  for  t(p,)  is  found  in  a  similar  manner.  Changing  vari¬ 
ables  in  (A1.7)  and  using  the  model  parameterization  in  (Al.9)  gives 

UTP 


r[Pi)  =2  J  (u2  -p,-2)1/2  (-^-)  du 


=  2  /  («2-p,?)i/2 


M 

Ew/  FAU) 


du 
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M 


-E-/ 


;'=* 


2/  («2  -  p,2)1/2  </« 


The  integral  is  evaluated  as 


(A1.13) 


2  f  (u2-  p?)1/2  du  =  u  (u2-??)1/2  -  p?  In[«  +  (^-p,2)1/2 


l 


J 


-  p,* 


ln[  tty  +(tty-^>^)1/2]  —  ln[  Uy_4  +  2)1/2  ] 


With  wmax  defined  ns  before  in  case  the  ray  turns  prior  to  the  interface. 
8.y  =  2  f  (u2-p?)1/2  du 


=  «max  (^ax-P.^)172  ~  «,_*  («/-,  -p,^)l/2 


Pi2 


«max  +(«max-P,'2)1/2 

«y-t  +  («/U-P,2)1/2 


(A1.14) 


Equation  Al.13  is  the  matrix  equation  for  r(p)  data,  7  =  Gw,  with  elements  of  G 
given  by  (Al.14). 

T(p)  may  be  found  by  re-arranging  the  top  equation  of  (Al.7). 

Tin)  =  PiAPi)+<Pi) 


Using  (Al.ll)  for  X(p,)  and  (Al.14)  for  7 ^p{)  results  in 


M 


T[Pi)  =  J]  «y 
2- i 


«max  («max-P.2)1/2  ~  «y_i  («/-!  -p2)1/2  +  Pi  In 


“max  +(«max-P,2)l/2 

My-!  +(«/_! -p,-2)1/2 


(A1.15) 

The  matrix  equation  is  T  =  Gw.  In  this  case  the  data  vector  consists  of  times 
T(pi),  and  gi}-  is  defined  by  the  quantity  in  the  outer  brackets. 
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In  each  case  (equations  Al.10,  Al.13,  and  Al.15)  the  medium  is  broken 
up  into  M  layers.  The  model  is  constructed  by  picking  values  of  slowness  (uy)  to 
specify  the  M+ 1  interfaces.  Slowness  at  the  surface  sets  the  first  interface 
(;'  =0).  The  last  interface  (/  =M)  should  be  set  slightly  below  the  turning  point 
of  the  deepest  ray  (i.e.  uu  should  be  less  than  the  smallest  ray  parameter  p  )  so 
that  the  ray  turns  in  the  last  layer. 

The  solution  contains  reciprocals  of  the  slowness  gradients  for  each  of 
the  M  layers.  Velocity  as  a  function  of  depth  is  recovered  by  first  finding  the 
layer  thicknesses  from  slowness  values  at  the  interfaces  and  the  slowness  gra¬ 
dients.  The  thickness  of  the  /**  layer  may  be  found  as 


(d*)i 


( du/dz)j 


=  («,■  -  «/_|)  Wy 


Depth  from  the  surface  to  any  interface  is  then  found  by  integrating  the  layer 
thicknesses  to  that  interface,  yielding  velocity  as  a  function  of  depth. 
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u(z)  = 


1 

v(z) 


depth  z 


Figure  A  1.1  Raypath  through  a  laterally  homogeneous  but  vertically  inhomogenous  medium. 

Slowness  is  assumed  to  decrease  with  depth  (or  equivalently,  velocity  increases 
with  depth).  Shows  the  ray’s  turning  point,  take-off  angle,  and  range.  Depth  at 
the  surface  is  z=0,  slowness  is  u(0). 


Figure  A1.2  Constant  gradient  model. 

Each  layer  has  a  constant  slowness  gradient. 
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Appendix  A2 

The  method  of  spectral  ratios 


Amplitude  loss  in  an  attenuating  medium 

Consider  a  signal  that  propagates  through  a  layer  of  an  attenuating 
medium  as  shown  in  Figure  A2.1.  The  pathlength  (layer  thickness  in  this  case) 
is  given  by  s.  Let  the  amplitude  spectrum  of  the  source  be  given  by  A0(f),  and 
the  amplitude  spectrum  of  the  received  time  series  be  A(f  ,s).  The  input  and 
output  are  related  by  the  transfer  function  H(f,s). 

Mf ,•)  -  A0(f)H(f  ,s)  (A2.i) 

H(f  ,s)  is  itself  made  up  of  a  number  of  functions  which  describe  the  changes  the 
signal  undergoes  as  it  propagates  through  the  layer.  In  general,  it  will  be  a 
function  of  the  pathlength  and  frequency.  The  received  waveform,  after  travel¬ 
ing  through  a  layer  of  sediment,  can  be  modeled  as 

Mf,*)  =  A0(f)H(f,s)  =  AM)[l(f)T{')Q(,)L(f,»)]  (A2.2) 


where 


A(/,s) 

is  the  amplitude  spectrum  of  the  received  signal 

Mf) 

is  the  amplitude  spectrum  of  the  source 

1(f) 

is. the  instrument  response 

represents  all  frequency  independent  losses 
(reflection  &  transmission  coefficients,  scattering,  etc.) 

<*(') 

represents  geometrical  spreading 

is  the  attenuation  operator 

(includes  all  frequency  dependent  attenuation) 

Attenuation  in  the  water  column  is  negligible  at  low  frequencies 
(pg.  106)],  so  a  good  representation  of  A0(f)  can  be  obtained  from  the 
water  path.  The  spectrum  of  the  direct  water  path  may  be  modeled  as 

[36,  65 
direct 

Ay)  (/  j  )  - 

(A2.3) 

where  <?„(«„)  is  geometrical  spreading  and  s„  is  the  pathlength  for  the  direct 
water  path.  By  dividing  the  spectrum  of  the  sediment  refracted  arrival  (A2.2) 
by  the  spectrum  of  the  direct  water  path  arrival  (A2.3),  factors  involving  the 
source  spectrum  and  instrument  response  are  eliminated. 
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A0(f) 


T 

s  H(f,s) 

k 


A(f,s) 


Figure  A2.1  Plane  wave  propagating  through  an  attenuating  medium. 


Z 


j-2 


Figure  A2.2  Ray  propagating  through  a  horizontally  layered,  attenuating  medium. 
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Systems  Command 
Washington.  D.  C.  20360 
Code  PD W- 124,  32QA 

Commander 

Naval  Ship  Res.  Sc  Dev.  Center 
Bethesda,  Maryland  20084 
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Code  L4Q,  L42 

Commanding  Officer 
Naval  Underwaeer  Systems  Center 
Newport,  Rhode  Island  20844 
Atten:  Dr.  K.  A.  Kemp,  Code  8213 

Officer  in  Charge 
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New  London  Laboratory 
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U.S.  Naval  Research  Laboratory 
Washington,  D,  C.  20375 
Code  2620,  2827,  5000,  5100,  5800 
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Silver  Spring,  Maryland  20910 
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University  of  Alaska 
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Director 
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Director 
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Department  of  Ocean  Engineering 
Florida  Atlantic  University- 
Boca  Raton,  Florida  33431 

Director 

Marine  Research  Laboratories 
c/o  Marine  Studies  Center 
University  of  Wisconsin 
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Director 
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P.O.  Box  30 

State  College,  Pennsylvania  16802 
Director 
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Director 
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Applied  Research  Laboratory 
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Director 
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Torrey  Cliff 

Palisades,  New  York  10964 
Director 

Woods  Hole  Oceanographic  Institution 
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Director 

Inst,  of  Ocean  Science  Engineering 
Catholic  University  of  America 
Washington.  D.C.  20017 

National  Science  Foundation 
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